[cig-commits] r14121 - seismo/3D/ADJOINT_TOMO/mtadj/utils

[liuqy at geodynamics.org]

Author: liuqy
Date: 2009-02-23 07:08:34 -0800 (Mon, 23 Feb 2009)
New Revision: 14121

Added:
   seismo/3D/ADJOINT_TOMO/mtadj/utils/combine_adj_src.pl
Modified:
   seismo/3D/ADJOINT_TOMO/mtadj/utils/readme
Log:
update readme file in utils directory
add combine_adj_src.pl script


Added: seismo/3D/ADJOINT_TOMO/mtadj/utils/combine_adj_src.pl
===================================================================
--- seismo/3D/ADJOINT_TOMO/mtadj/utils/combine_adj_src.pl	                        (rev 0)
+++ seismo/3D/ADJOINT_TOMO/mtadj/utils/combine_adj_src.pl	2009-02-23 15:08:34 UTC (rev 14121)
@@ -0,0 +1,90 @@
+#!/usr/bin/perl -w
+
+#--------------------------------------------------------------------
+# combine_adj_src.pl
+# Qinya Liu, 05-May-2007
+#
+# This script combines two sets of adjoint sources (for example, P wave time windows
+# and surface wave time windows), and outputs the new adjoint sources, along with
+# a new STATIONS file to be used.
+#
+# To remove the cc label on files:
+#   mv_files.pl "*.cc.adj" "*.adj"
+#
+# CALLED BY: combine_adj_src_all.pl
+#
+#--------------------------------------------------------------------
+
+if (@ARGV != 5) {die("Usage: combine_adj_src.pl DIR_1 DIR_2 DIR_NEW type1(cc,ik,mtm) type2(cc,ik,mtm)\n");}
+
+$dir1 = $ARGV[0];
+$dir2 = $ARGV[1];
+
+$dir_new = $ARGV[2];
+$type1 = $ARGV[3];
+$type2 = $ARGV[4];
+
+$sta1 = "$dir1/STATIONS_ADJOINT";
+$sta2 = "$dir2/STATIONS_ADJOINT";
+$sta_new = "$dir_new/STATIONS_ADJOINT";
+if (not -f $sta1 or not -f $sta2) {die("Check if $sta1 and $sta2 exist or not\n");}
+if (not -d $dir_new) {system("mkdir $dir_new");}
+
+# open STATIONS file 1
+open(FILE,"$sta1");
+(@sta) = <FILE>;    # first line is the number of stations
+close(FILE);
+for ($i=1; $i<@sta; $i++){
+  ($sta,$nt,$lat,$lon,$dep,$bur) = split(" ",$sta[$i]);
+  $adj{$sta}{num} = 1;
+  $adj{$sta}{net} = "$nt";
+  $adj{$sta}{info} = "$nt $lat $lon $dep $bur";
+  @files = glob("$dir1/$sta.$nt.BH?.*adj");
+  if (@files == 0) {print("Check if $dir1/$sta.$nt.BH?.*adj exist or not\n");}
+}
+
+# open STATIONS file 2
+open(FILE,"$sta2");
+(@sta) = <FILE>;    # first line is the number of stations
+close(FILE);
+for ($i=1; $i<@sta; $i++){
+  ($sta,$nt,$lat,$lon,$dep,$bur) = split(" ",$sta[$i]);
+  if (not defined $adj{$sta}{num}) {
+    $adj{$sta}{num} = 1;
+    $adj{$sta}{net} = "$nt";
+    $adj{$sta}{info} = "$nt $lat $lon $dep $bur";}
+  else {$adj{$sta}{num} ++ ; 
+        # for example, this crashes if HAST.TA and HAST.BK are included
+	if ($nt ne $adj{$sta}{net}) {die("Check if network name same for $sta\n");}
+  }
+  @files = glob("$dir2/$sta.$nt.BH?.*adj");
+  if (@files == 0) {print("Check if $dir2/$sta.$nt.BH?.*adj exist or not\n");}
+}
+
+
+open(FILE,">$dir_new/sta.tmp"); $nsta = 0;
+foreach $sta (keys %adj) {
+#foreach $sta ("TCF") {
+  $net = $adj{$sta}{net}; $nsta ++ ;
+  print "$nsta, $sta, $net\n";
+  print FILE "$sta $adj{$sta}{info} \n";
+  if ($adj{$sta}{num} > 1) {
+    # add up two adjoint sources
+    foreach $comp ("BHE","BHN","BHZ") {
+      system("paste $dir1/$sta.$net.$comp.$type1.adj $dir2/$sta.$net.$comp.$type2.adj | awk '{print \$1, \$2+\$4}' > $dir_new/$sta.$net.$comp.adj");
+#      print "paste $dir1/$sta.$net.$comp.$type1.adj $dir2/$sta.$net.$comp.$type2.adj | awk '{print \$1, \$2+\$4}' > $dir_new/$sta.$net.$comp.adj\n";
+    }
+  }else{
+    # copy over the adjoint source
+    foreach $comp ("BHE","BHN","BHZ") {
+       system("cp -f $dir1/$sta.$net.$comp.$type1.adj  $dir_new/$sta.$net.$comp.adj >& /dev/null");
+       system("cp -f $dir2/$sta.$net.$comp.$type2.adj  $dir_new/$sta.$net.$comp.adj >& /dev/null");
+     }
+  }
+}
+close(FILE);
+
+print "\n number of stations in new STATIONS file should be $nsta\n";
+system(" echo $nsta > $sta_new; cat $dir_new/sta.tmp >> $sta_new");
+
+#----------------------------


Property changes on: seismo/3D/ADJOINT_TOMO/mtadj/utils/combine_adj_src.pl
___________________________________________________________________
Name: svn:executable
   + *

Modified: seismo/3D/ADJOINT_TOMO/mtadj/utils/readme
===================================================================
--- seismo/3D/ADJOINT_TOMO/mtadj/utils/readme	2009-02-23 07:39:57 UTC (rev 14120)
+++ seismo/3D/ADJOINT_TOMO/mtadj/utils/readme	2009-02-23 15:08:34 UTC (rev 14121)
@@ -1,3 +1,35 @@
-rotate_adj_src utility helps to rotate the adjoint sources
+1. rotate_adj_src utility helps to rotate the adjoint sources
 calculated from the measurements to E, N, Z components which
 is ready to go directly into SEM adjoint code.
+
+2. combine_adj_src.pl -- directly copied over from old version. I assume
+     the original purpose is to combine cc adjoint source for body
+     waves and mt adjoint source for surface waves, it is probably
+     not the best thing to use in practice, considering the possible
+     weight one needs to add to the body waves in actual inversions.
+
+2. A list of scripts that might worth transfering from old version
+(a lot of them by Carl), notably run_mt_measure_adj.pl and 
+plot_win_adj.pl. Keep in mind now our output seismograms are in SAC format.
+
+./scripts_meas/mtm.pl
+./scripts_meas/plot_adj.pl
+./scripts_meas/plot_mtm.pl
+./scripts_meas/plot_cc.pl
+
+./scripts_tomo/combine_2_adj_src.pl
+./scripts_tomo/run_tomo.pl
+./scripts_tomo/prepare_mt_measure_adj.pl
+./scripts_tomo/combine_3_adj_src_all.pl
+./scripts_tomo/combine_3_adj_src.pl
+./scripts_tomo/run_mt_cc_plot.pl
+./scripts_tomo/combine_2_adj_src_all.pl
+./scripts_tomo/run_mt_measure_adj.pl
+
+./PLOTS/make_pdf_by_event.pl
+./PLOTS/make_pdf_by_station.pl
+./PLOTS/plot_win_adj_all.pl
+./PLOTS/plot_win_stats_all.pl
+./PLOTS/plot_win_stats.pl
+./PLOTS/plot_win_adj.pl
+