[cig-commits] r13643 - seismo/3D/GRD_CMT3D/scripts
liuqy at geodynamics.org
liuqy at geodynamics.org
Wed Dec 10 13:07:27 PST 2008
Author: liuqy
Date: 2008-12-10 13:07:27 -0800 (Wed, 10 Dec 2008)
New Revision: 13643
Added:
seismo/3D/GRD_CMT3D/scripts/cp_syn.bash
seismo/3D/GRD_CMT3D/scripts/mech_table2.pl
seismo/3D/GRD_CMT3D/scripts/plot_socal_events.pl
Removed:
seismo/3D/GRD_CMT3D/scripts/cp_syn_all.bash
Modified:
seismo/3D/GRD_CMT3D/scripts/mech_table.pl
seismo/3D/GRD_CMT3D/scripts/plot_data_and_syn.pl
seismo/3D/GRD_CMT3D/scripts/plot_grid.pl
seismo/3D/GRD_CMT3D/scripts/process_data_and_syn.pl
Log:
Update scripts:
plot_data_and_syn.pl -- add capability of plotting new synthetics
mech_table2.pl -- inversion table for 2,3,6-30 sec
plot_socal_events.pl -- socal topo plot
Added: seismo/3D/GRD_CMT3D/scripts/cp_syn.bash
===================================================================
--- seismo/3D/GRD_CMT3D/scripts/cp_syn.bash (rev 0)
+++ seismo/3D/GRD_CMT3D/scripts/cp_syn.bash 2008-12-10 21:07:27 UTC (rev 13643)
@@ -0,0 +1,72 @@
+#!/bin/bash
+
+if [ $# != 1 ]; then
+ echo "Usage: cp_syn.bash evid-file"; exit
+fi
+
+extra_dir=""
+
+for evid in `cat $1 `; do
+ echo "==== $evid ===="
+ mkdir -p ${evid}${extra_dir};
+ cd ${evid}${extra_dir}; mkdir -p syn; mkdir -p backup
+## copy sem syn from cluster
+# scp $cluster/${evid}${extra_dir}/* syn/
+ mv ../syns/*${evid}*.bz syn/
+ cp ../event_info/CMTSOLUTION_$evid CMTSOLUTION; cp CMTSOLUTION syn/
+ cp ../event_info/STATIONS .
+## Unzip synthetics
+ cd syn
+ tar xjvf CMTSOLUTION_${evid}.tar.bz
+ for ext in Mrr Mtt Mpp Mrt Mrp Mtp dep lat lon; do
+ echo " Untar $ext ..."; tar xjvf ${evid}_${ext}.tar.bz >/dev/null
+ if [ $? != 0 ]; then
+ echo "Error untaring ... $ext ..."; exit
+ fi
+ mv -f ${evid}_${ext}.tar.bz ../backup
+ rename "s/semd.sac.$ext/$ext/" *.semd.sac.$ext
+ done
+ ls -1 *.Mrr | perl -pi -e 's/\.Mrr//g' -- > file_list
+ xadd_frechet_derivatives s file_list ../CMTSOLUTION 1.0e22
+ cd ..
+## Copy measurement file
+ cp ../event_info/MEASUREMENT_WINDOWS_${evid}_T00?_T030_m12 .
+## here check if MEASURE file has the right number of data/syn pairsa
+ mafile="MEASUREMENT_WINDOWS_${evid}_ALL"
+ nma=0
+ rm -f $mafile
+ for freq in 2 3 6; do
+ mfile="MEASUREMENT_WINDOWS_${evid}_T00${freq}_T030_m12"
+ n0=`grep "^[DS]" $mfile | wc | awk '{print $1}'`
+ n1=`echo "$n0 / 2" | bc`
+ n2=`head -1 $mfile`
+ if [ $n1 != $n2 ]; then
+ echo "Number of files inconsistent in $mfile: $n1 and $n2"; exit
+ fi
+## here modify dir names: DATA -> data_T006_T030, SYN -> syn_T006_T030
+ perl -pi -e "s/DATA/data_T00${freq}_T030/g" $mfile
+ perl -pi -e "s/SYN/syn_T00${freq}_T030/g" $mfile
+ perl -pi -e "s/\.T00${freq}_T030//g" $mfile
+ perl -pi -e "s/\.semd\.sac\.m12//g" $mfile
+ awk 'NR > 1 {print $0}' $mfile >> $mafile
+ nm0=`head -n 1 $mfile`
+ nma=`echo $nm0 + $nma | bc `
+ done
+## right number of measurements
+ echo $nma | cat - $mafile > out.tmp; mv out.tmp $mafile
+#### Link cmt3d and grid3d_flexwin
+ ln -s ../../grd_cmt3d/cmt3d/cmt3d_flexwin
+ ln -s ../../grd_cmt3d/grid3d/grid3d_flexwin
+## here modify INVERSION.PAR GRID3D.PAR for MEASURE file
+ cp ../../grd_cmt3d/cmt3d/INVERSION.PAR .
+ cp ../../grd_cmt3d/grid3d/GRID3D.PAR .
+ perl -pi -e "s/flexwin\.out/$mafile/g" INVERSION.PAR
+ cp INVERSION.PAR INVERSION.PAR.SAVE
+ perl -pi -e "s/flexwin\.out/$mafile/g" GRID3D.PAR
+## make sure you have the right dmoment, ddepth, and dlocation
+ dlat=`diff -b syn/CMTSOLUTION syn/CMTSOLUTION_lat | grep lat | awk '{print $3}' | perl -pi -e 's/\n/ /g' | awk '{print $2-$1}'`
+ ddep=`diff -b syn/CMTSOLUTION syn/CMTSOLUTION_dep | grep dep | awk '{print $3}' | perl -pi -e 's/\n/ /g' | awk '{print $2-$1}'`
+ perl -pi -e "s/^.*1.0e22*/$dlat $ddep 1.0e22/" INVERSION.PAR
+####
+ cd ..;
+done
Property changes on: seismo/3D/GRD_CMT3D/scripts/cp_syn.bash
___________________________________________________________________
Name: svn:executable
+ *
Deleted: seismo/3D/GRD_CMT3D/scripts/cp_syn_all.bash
===================================================================
--- seismo/3D/GRD_CMT3D/scripts/cp_syn_all.bash 2008-12-10 13:50:34 UTC (rev 13642)
+++ seismo/3D/GRD_CMT3D/scripts/cp_syn_all.bash 2008-12-10 21:07:27 UTC (rev 13643)
@@ -1,73 +0,0 @@
-#!/bin/bash
-
-if [ $# != 1 ]; then
- echo "Usage: cp_syn.bash evid-file"; exit
-fi
-
-cluster="lqy@$wag:/net/wagholi/scratch1/lqy/CMT3D_Basin_Adjoint/cmt_mass/"
-
-extra_dir="_5.6km"
-
-for evid in `cat $1 `; do
- mkdir -p ${evid}${extra_dir};
-## copy sem syn from cluster
- cd ${evid}${extra_dir}; mkdir -p syn; mkdir -p backup
- scp $cluster/${evid}${extra_dir}/* syn/
- cp syn/CMTSOLUTION .; cp ../event_info/STATIONS .
-
- cp ../event_info/MEASUREMENT_WINDOWS_${evid}_T00?_T030_m12 .
-## here check if MEASURE file has the right number of data/syn pairsa
- mafile="MEASUREMENT_WINDOWS_${evid}_ALL"
- nma=0
- rm -f $mafile
- for freq in 2 3 6; do
- mfile="MEASUREMENT_WINDOWS_${evid}_T00${freq}_T030_m12"
- n0=`grep "^[DS]" $mfile | wc | awk '{print $1}'`
- n1=`echo "$n0 / 2" | bc`
- n2=`head -1 $mfile`
- if [ $n1 != $n2 ]; then
- echo "Number of files inconsistent in $mfile: $n1 and $n2"; exit
- fi
-## here to modify DATA -> data_T006_T030, SYN -> syn_T006_T030
- perl -pi -e "s/DATA/data_T00${freq}_T030/g" $mfile
- perl -pi -e "s/SYN/syn_T00${freq}_T030/g" $mfile
- perl -pi -e "s/\.T00${freq}_T030//g" $mfile
- perl -pi -e "s/\.semd\.sac\.m12//g" $mfile
- awk 'NR > 1 {print $0}' $mfile >> $mafile
- nm0=`head -n 1 $mfile`
- nma=`echo $nm0 + $nma | bc `
- rm -f $mfile
- done
-## here to link cmt3d/grid3d_flexwin
- ln -s ../../grd_cmt3d/cmt3d/cmt3d_flexwin
- ln -s ../../grd_cmt3d/grid3d/grid3d_flexwin
-## here to modify INVERSION.PAR GRID3D.PAR for MEASURE file
- cp ../../grd_cmt3d/cmt3d/INVERSION.PAR .
- cp ../../grd_cmt3d/grid3d/GRID3D.PAR .
- perl -pi -e "s/flexwin\.out/$mafile/g" INVERSION.PAR
- cp INVERSION.PAR INVERSION.PAR.SAVE
- perl -pi -e "s/flexwin\.out/$mafile/g" GRID3D.PAR
-## make sure you have the right dmoment, ddepth, and dlocation
- dlat=`diff -b syn/CMTSOLUTION syn/CMTSOLUTION_lat | grep lat | awk '{print $3}' | perl -pi -e 's/\n/ /g' | awk '{print $2-$1}'`
- ddep=`diff -b syn/CMTSOLUTION syn/CMTSOLUTION_dep | grep dep | awk '{print $3}' | perl -pi -e 's/\n/ /g' | awk '{print $2-$1}'`
- perl -pi -e "s/^.*1.0e22*/$dlat $ddep 1.0e22/" INVERSION.PAR
- perl -pi -e "s/^.*1.0e22*/$dlat $ddep 1.0e22/" GRID3D.PAR
-
-## right number of measurements
- echo $nma | cat - $mafile > out.tmp; mv out.tmp $mafile
-## Unzip synthetics
- cd syn
- for ext in Mrr Mtt Mpp Mrt Mrp Mtp dep lat lon; do
- tar xjvf ${evid}_${ext}.tar.bz >/dev/null
- if [ $? != 0 ]; then
- echo "Error untaring ... $ext ..."; exit
- fi
- mv -f ${evid}_${ext}.tar.bz ../backup
- rename "s/semd.sac.$ext/$ext/" *.semd.sac.$ext
- done
- ls -1 *.Mrr | perl -pi -e 's/\.Mrr//g' -- > file_list
- xadd_frechet_derivatives s file_list ../CMTSOLUTION 1.0e22
- cd ..
-####
- cd ..;
-done
Modified: seismo/3D/GRD_CMT3D/scripts/mech_table.pl
===================================================================
--- seismo/3D/GRD_CMT3D/scripts/mech_table.pl 2008-12-10 13:50:34 UTC (rev 13642)
+++ seismo/3D/GRD_CMT3D/scripts/mech_table.pl 2008-12-10 21:07:27 UTC (rev 13643)
@@ -59,6 +59,7 @@
print BASH "cp -f INVERSION.$k INVERSION.PAR\n";
print BASH "cmt3d_flexwin > cmt3d.stdout\n";
print BASH "if [ \$? != 0 ]; then\n echo \"exiting case $k\";exit\nfi\n";
+ print BASH "perl -pi -e \"s/INV/$names[$k]/\" $cmt_new\n";
print BASH "mv -f $cmt_new $outdir/$cmt[$k]\n";
print BASH "mv -f cmt3d_flexwin.out $outdir/cmt3d_out.$k\n";
print BASH "mv -f cmt3d.stdout $outdir/cmt3d_stdout.$k\n";
Added: seismo/3D/GRD_CMT3D/scripts/mech_table2.pl
===================================================================
--- seismo/3D/GRD_CMT3D/scripts/mech_table2.pl (rev 0)
+++ seismo/3D/GRD_CMT3D/scripts/mech_table2.pl 2008-12-10 21:07:27 UTC (rev 13643)
@@ -0,0 +1,123 @@
+#!/usr/bin/perl
+use lib '/opt/seismo/lib/perl';
+use GMT_PLOT;
+use GMT_PLOT_SC;
+use CMT_TOOLS;
+use POSIX;
+
+## ok, this version includes both 2, 3 and 6 seconds measurements
+# parameter files
+$inv_file="INVERSION.PAR.SAVE";
+$outdir="mech"; $psfile="mech.ps";
+
+$weigh_data="2 2 1 0.5 1.15 0.55 0.78";
+$weigh_data_nocomp="1 1 1 0.5 1.15 0.55 0.78";
+$weigh_data_noaz="2 2 1 0 1.15 0.55 0.78";
+
+if (not -f "CMTSOLUTION") {die("Check cmt solution file\n");}
+($evid)=split(" ",`cat CMTSOLUTION | awk 'NR == 2 {print \$3}'`);
+$mwin="MEASUREMENT_WINDOWS";
+
+ at names=("Initial","Grid Search","7 Par+ZT",
+ "7 Par+no-SC","6 par+DC","7 Par+DC",
+ "7 Par(2-30s)","7 Par(3-30s)","7 Par(6-30s)");
+ at npars=(6,6,7, 7,6,7, 7,7,7);
+ at wdatas=($weigh_data, $weigh_data, $weigh_data, $weigh_data,$weigh_data,$weigh_data, $weigh_data,$weigh_data,$weigh_data);
+ at scorrs=("true","true",".true.", ".false.", ".true.", ".true.", ".true.",".true.",".true.");
+ at cons=(".true. true. 0.0", ".true. .true. 0.0", ".true. .false. 0.0",
+ ".true. .false. 0.0", "true. .true. 0.0", "true. .true. 0.0",
+ ".true. .true. 0.0",".true. .true. 0.0",".true. .true. 0.0");
+ at meas=("${mwin}_${evid}_ALL","${mwin}_${evid}_ALL","${mwin}_${evid}_ALL",
+ "${mwin}_${evid}_ALL","${mwin}_${evid}_ALL","${mwin}_${evid}_ALL",
+ "${mwin}_${evid}_T002_T030_m12","${mwin}_${evid}_T003_T030_m12","${mwin}_${evid}_T006_T030_m12");
+
+$nbeach=@npars;
+$ncols = 3;
+$nrows = 3;
+# plot size jx/jy, plot origin in xy
+($jx,$jy) = shift_xy("1 1","$ncols $nrows","7 8", \@xy,"0.8 0.8");
+$JX = "-JX$jx/$jy"; $R="-R0/3/0/3";
+$GMT_PLOT::paper_orient="-P";
+
+if (not -f "cmt3d_flexwin") {die("Check if cmt3d_flexwin exists or not\n");}
+
+open(INV,"$inv_file")|| die("check $inv_file\n");
+ at all=<INV>;close(INV);
+
+$cmt=$all[0]; $cmt_new=$all[1]; $npar=$all[2]; $delta=$all[3];
+$flex=$all[4];$lwdata=$all[5]; $wdata=$all[6]; $scorr=$all[7]; $con=$all[8];
+$wns=$all[9];
+chomp($cmt_new);chomp($cmt);
+(undef,undef,$ename)=split(" ",`grep 'event name' $cmt`);
+
+ at cmt=($cmt);@cmt=(@cmt,"${cmt}_GRD");
+if (-f "${cmt}_GRD") {print BASH "cp -f ${cmt}_GRD $outdir\n";}
+for ($i=2;$i<$nbeach;$i++) {@cmt=(@cmt,"CMTSOLUTION.$i");}
+
+open(BASH,">mech_table.bash");
+print BASH "mkdir -p $outdir\n";
+print BASH "cp -f $cmt $outdir\n";
+if (-f "${cmt}_GRD") {print BASH "cp -f ${cmt}_GRD $outdir\n";}
+
+for ($k=2;$k<$nbeach;$k++) { # cmt3d runs
+ print BASH "echo Running cmt3d_flexwin for INVERSION.$k ...\n";
+ open(INV,">INVERSION.$k");
+ print INV "${cmt}\n${cmt_new}\n$npars[$k]\n$delta$meas[$k]\n.true.\n$wdatas[$k]\n$scorrs[$k]\n$cons[$k]\n";
+ if ($k==2) {print INV ".true.\n";} else { print INV "${wns}\n";}
+ close(INV);
+ print BASH "# --- mech $k ----\n";
+ print BASH "cp -f INVERSION.$k INVERSION.PAR\n";
+ print BASH "cmt3d_flexwin > cmt3d.stdout\n";
+ print BASH "if [ \$? != 0 ]; then\n echo \"exiting case $k\";exit\nfi\n";
+ print BASH "perl -pi -e \"s/INV/$names[$k]/\" $cmt_new\n";
+ print BASH "mv -f $cmt_new $outdir/$cmt[$k]\n";
+ print BASH "mv -f cmt3d_flexwin.out $outdir/cmt3d_out.$k\n";
+ print BASH "mv -f cmt3d.stdout $outdir/cmt3d_stdout.$k\n";
+ print BASH "mv -f INVERSION.$k $outdir\n";
+}
+close(BASH);
+system("chmod a+x mech_table.bash; mech_table.bash");
+for ($k=2;$k<$nbeach;$k++) {
+ (@tmp)=split(" ",`grep Variance $outdir/cmt3d_stdout.$k`);
+ $var[$k]=$tmp[8]; print "$k -- VR = $var[$k]\n";}
+#die("here\n");
+
+
+open(BASH,">mech_plot.bash");
+print BASH "gmtset BASEMAP_TYPE plain ANOT_FONT_SIZE 9 HEADER_FONT_SIZE 10 MEASURE_UNIT inch PAPER_MEDIA letter TICK_LENGTH 0.1c\n";
+
+plot_psxy(\*BASH,$psfile,"$JX -K -X0 -Y0 $R","");
+$k=0;
+for ($i=0;$i<$nrows;$i++) {
+ for ($j=0;$j<$ncols;$j++) {
+ print BASH "# ---- Mech $k, $cmt[$k]------\n";
+ print "---- Mech $k, $cmt[$k]------\n";
+ $xy=$xy[$i][$j];($x,$y) = split(/\//,$xy);
+ plot_psxy(\*BASH,$psfile,"$JX -X$x -Y$y $R","");
+ if (-f "$outdir/$cmt[$k]") {
+ plot_psmeca_raw(\*BASH,$psfile,"-JX -R -Sm1.0 -B3/3wesn","1.5 2","$outdir/$cmt[$k]");
+ @output = `cmtsol2faultpar.pl $outdir/$cmt[$k]`;
+ (undef,$tmp1) = split(" ", $output[4]);
+ (undef,$tmp2) = split(" ", $output[6]);
+ ($s,$d,$r,$p,$t) = split(/\//,$tmp2);
+ $tex="1 1 9 0 4 CM $names[$k]\n";
+ $tex.="0.2 0.8 9 0 4 LM Mw/dep/eps/= $tmp1\n0.2 0.6 9 0 4 LM S/D/R= $s/$d/$r";
+ if ($k>=2) {$tex.=sprintf("\n0.2 0.4 9 0 4 LM Var=%6.2f%",$var[$k]);}
+ plot_pstext(\*BASH,$psfile,"-JX -R -B -G0/0/255","$tex");
+ }
+ if ($k==1) {
+ ($nt) = split(" ",`grep 'T.sac.d' $meas[2] | wc | awk '{print \$1}'`);
+ ($nr) = split(" ",`grep 'R.sac.d' $meas[2] | wc | awk '{print \$1}'`);
+ ($nz) = split(" ",`grep 'Z.sac.d' $meas[2] | wc | awk '{print \$1}'`);
+ ($n2) = split(" ",`grep 'T002.*.sac.d' $meas[2] | wc | awk '{print \$1}'`);
+ ($n3) = split(" ",`grep 'T003.*.sac.d' $meas[2] | wc | awk '{print \$1}'`);
+ ($n6) = split(" ",`grep 'T006.*.sac.d' $meas[2] | wc | awk '{print \$1}'`);
+
+ plot_pstext(\*BASH,$psfile,"-JX -R -B -N","1.5 4 12 0 4 CM $ename ($nz/$nr/$nt; $n2/$n3/$n6)");}
+ plot_psxy(\*BASH,$psfile,"-JX -X-$x -Y-$y","");
+ $k++;
+ }
+}
+plot_psxy(\*BASH,$psfile,"-JX -O","");
+close(BASH);
+system("chmod a+x mech_plot.bash; mech_plot.bash");
Property changes on: seismo/3D/GRD_CMT3D/scripts/mech_table2.pl
___________________________________________________________________
Name: svn:executable
+ *
Modified: seismo/3D/GRD_CMT3D/scripts/plot_data_and_syn.pl
===================================================================
--- seismo/3D/GRD_CMT3D/scripts/plot_data_and_syn.pl 2008-12-10 13:50:34 UTC (rev 13642)
+++ seismo/3D/GRD_CMT3D/scripts/plot_data_and_syn.pl 2008-12-10 21:07:27 UTC (rev 13643)
@@ -14,9 +14,9 @@
use List::Util qw[min max];
-if (@ARGV == 0) {die("Usage: plot_data_and_syn.pl -d data_dir,data_ext -s syn_dir,syn_ext -M flexwin_file -m CMTSOLUTION \n");}
+if (@ARGV == 0) {die("Usage: plot_data_and_syn.pl -d data_dir,data_ext -s syn_dir,syn_ext -Snew_syn_dir,new_syn_ext -M flexwin_file -m CMTSOLUTION -A Z/R/T(0.03/0.03/0.01)\n");}
-if (!getopts('d:s:m:M:')) {die("Check input options\n");}
+if (!getopts('d:s:S:m:M:A:')) {die("Check input options\n");}
# data/syn dir and extension
if ($opt_M) {
@@ -33,8 +33,13 @@
else {$syn_ext="";}}
} else {die("Input -M or -d/-s options\n");}
+if ($opt_S) {$plot_new = 1;
+ ($new_syn_dir,$new_syn_ext) = split(/\,/,$opt_S);
+ if ($new_syn_ext) {$new_syn_ext=".$new_syn_ext";}
+ else {$new_syn_ext="";}}
+
if ($opt_m) {$cmt_file=$opt_m;} else {die("Input cmt file \n");}
-($moment) = split(" ",` cmtsol2faultpar.pl $cmt_file | awk 'NR == 3 {print \$3}'`);
+($moment) = split(" ",`cmtsol2faultpar.pl $cmt_file | awk 'NR == 3 {print \$3}'`);
# comps
$ncols=3; @comps=("Z","R","T");#Z R T
@@ -52,11 +57,16 @@
$Vl=3.0; # surface wave velocity to determine tend
$bp=10; $bp2=10;# start time before t1 header and end time after tend
$E="-Ent1"; # align by trace number and t1
- at size_all=(0.03/2.4e23*$moment,0.03/2.4e23*$moment,0.06/2.4e23*$moment);
+if ($opt_A) {
+ ($za,$ra,$ta) = split(/\//,$opt_A);
+ @size_all=($za/2.4e23*$moment,$ra/2.4e23*$moment,$ta/2.4e23*$moment);
+} else {
+ @size_all=(0.03/2.4e23*$moment,0.03/2.4e23*$moment,0.06/2.4e23*$moment);}
print "plot sizes = @size_all \n";
-$wdata="-W3/0/0/0"; $wsyn="-W2/255/0/0"; $gwin="-G220 -W3";
+$wdata="-W3/0/0/0"; $wsyn="-W2/255/0/0"; $gwin="-G220 -W3";
+if ($plot_new) {$wsyn="-W2/0/0/255"; $wnsyn="-W2/255/0/0";}
- $jxy_map="-X3 -Y8"; $JM="-JM3";
+$jxy_map="-X3 -Y8"; $JM="-JM3";
# *****************************************************
# loop over data files, select corresponding syn files and set up % file
@@ -68,16 +78,19 @@
(undef,$net,$sta,$chan,$dist,$az,$t1,$stla,$stlo) = split(" ",`saclst knetwk kstnm kcmpnm dist az t1 stla stlo f $datfile`);
($comp) = split(" ",`echo $chan | awk '{print substr(\$1,3,1)}'`);
$synfile = "${syn_dir}/${sta}.${net}.BH${comp}${syn_ext}";
+
if (not -f $synfile) {
print "**** Missing synthetic pair $synfile for $datfile ***, skipping ...\n"; }
else {
$file{$sta}{dist}=$dist; $file{$sta}{az}=$az; $file{$sta}{P}=$t1;
$file{$sta}{"data-$comp"}=$datfile; $file{$sta}{"syn-$comp"}=$synfile;
$file{$sta}{stla}=$stla; $file{$sta}{stlo}=$stlo;
+ if ($plot_new) {
+ $new_synfile="${new_syn_dir}/${sta}.${net}.BH${comp}${new_syn_ext}";
+ $file{$sta}{"new-syn-$comp"}=$new_synfile;}
}
}
} else { # given flexwin output
-# (@files) = `grep -v '^ *-*[0-9].*$' $mfile`;
open(MFILE,"$mfile");
$nfile=<MFILE>;
for ($i=0;$i<$nfile;$i++) {
@@ -92,6 +105,9 @@
$file{$sta}{dist}=$dist; $file{$sta}{az}=$az; $file{$sta}{P}=$t1;
$file{$sta}{"data-$comp"}=$datfile; $file{$sta}{"syn-$comp"}=$synfile;
$file{$sta}{stla}=$stla; $file{$sta}{stlo}=$stlo;
+ if ($plot_new) {
+ $new_synfile="${new_syn_dir}/${sta}.${net}.BH${comp}${new_syn_ext}";
+ $file{$sta}{"new-syn-$comp"}=$new_synfile;}
# windows
($nwin)=split(" ",<MFILE>);
$file{$sta}{"nwin-$comp"}=$nwin; $file{$sta}{"win-$comp"}="";
@@ -146,11 +162,12 @@
# ***** CONTROL PAR ******
$size=2*$size_all[$i]/($dist_min+$dist_max);
- $data_files="";$syn_files="";$win_files="";
+ $data_files="";$syn_files="";$win_files=""; $new_syn_files="";
if ($i==0) {$tex_file="";$xy_file="";$t_file="";}
for ($k=0;$k<@stas;$k++) {
$dfile=$file{$stas[$k]}{"data-$comps[$i]"};
$sfile=$file{$stas[$k]}{"syn-$comps[$i]"};
+ if ($plot_new) {$nsfile=$file{$stas[$k]}{"new-syn-$comps[$i]"};}
# station name, dist, az
if ($i==0) {
$kt1=$k+0.2; $tt3=$tmin; $kt2 = $k-0.2;
@@ -165,6 +182,7 @@
if (defined $dfile and defined $sfile) {
$data_files.="$dfile 0 $k\n"; $syn_files.="$sfile 0 $k\n";
+ if ($plot_new) {$new_syn_files.="$nsfile 0 $k\n";}
# windows
if ($flexwin) {
@wins=split("\n",$file{$stas[$k]}{"win-$comps[$i]"});
@@ -176,7 +194,7 @@
}
}
}
- chomp($data_files); chomp($syn_files);
+ chomp($data_files); chomp($syn_files); chomp($new_syn_files);
# plot windows
if ($flexwin) {
@@ -186,7 +204,10 @@
# plot data/syn
plot_pssac2_raw(\*BASH,$psfile,"$JX $R $B[$i] $E $C -M$size/0 $wdata -N","$data_files");
plot_pssac2_raw(\*BASH,$psfile,"-JX -R -B $E $C -M$size/0 $wsyn -N","$syn_files");
- if ($i==0) {
+ if ($plot_new) {
+ plot_pssac2_raw(\*BASH,$psfile,"-JX -R -B $E $C -M$size/0 $wnsyn -N","$new_syn_files");
+ }
+ if ($i==0) {
chomp($tex_file); chomp($t_file); chomp($xy_file);
plot_pstext(\*BASH,$psfile,"-JX -R -B","$tex_file");}
Modified: seismo/3D/GRD_CMT3D/scripts/plot_grid.pl
===================================================================
--- seismo/3D/GRD_CMT3D/scripts/plot_grid.pl 2008-12-10 13:50:34 UTC (rev 13642)
+++ seismo/3D/GRD_CMT3D/scripts/plot_grid.pl 2008-12-10 21:07:27 UTC (rev 13643)
@@ -83,7 +83,7 @@
$arr=$ar[-1]+$ddr*2; # position for title
$pp="$as[$nrows/2] $arr";
- print CSH "makecpt -T$min/$max/$db2 -Cseis -Z -I > temp.cpt \n";
+ print CSH "makecpt -T$min/$max/$db2 -Cseis -Z > temp.cpt \n";
plot_psxy(\*CSH,$psfile,"$JX -K -X0 -Y0","");
Added: seismo/3D/GRD_CMT3D/scripts/plot_socal_events.pl
===================================================================
--- seismo/3D/GRD_CMT3D/scripts/plot_socal_events.pl (rev 0)
+++ seismo/3D/GRD_CMT3D/scripts/plot_socal_events.pl 2008-12-10 21:07:27 UTC (rev 13643)
@@ -0,0 +1,44 @@
+#!/usr/bin/perl
+
+# this script plots given station and cmt files on
+# a southern california topography map.
+# Special thanks to Carl's master script and his effort
+# of putting together the topo-bathymetry grid file for
+# this region
+
+use lib '/opt/seismo/lib/perl';
+use GMT_PLOT;
+use GMT_PLOT_SC;
+use Getopt::Std;
+
+if (@ARGV == 0) {die("plot_socal_events.pl sta-file cmt-solution-files\n");}
+
+$station=$ARGV[0];
+if (not -f $station) {die("Check $station file\n");}
+ at cmt_files=@ARGV[1..$#ARGV];
+foreach $file (@cmt_files) {
+ if (not -f $file) {die("Check cmt file $file\n");}
+}
+
+# this controls the paper orientation;
+
+$file="/data2/Datalib/SC/w140n40.Bathmetry.srtm.swap";
+$R="-122/-114.5/32/37";
+$JM="-JM7.5i";
+$psfile="socal_map.ps";
+
+open(BASH,">socal_map.bash");
+
+print BASH "gmtset PAPER_MEDIA letter BASEMAP_TYPE plain PLOT_DEGREE_FORMAT D TICK_LENGTH 0.3c LABEL_FONT_SIZE 12 ANOT_FONT_SIZE 10 HEADER_FONT 1 ANOT_FONT 1 LABEL_FONT 1 HEADER_FONT_SIZE 18 FRAME_PEN 2p TICK_PEN 2p \n";
+
+plot_sc_all(\*BASH,$psfile,"$JM -R$R -W1.5 -w0.5 -B2/2WesN -W1p,0 -Na/1.0p,255/255/255 -I$file.int -C$file.cpt $file.grd -K","");
+
+plot_psxy(\*BASH,$psfile,"-JM -R -St0.07 -W0.3 -G255/0/0","$station");
+
+plot_psmeca(\*BASH,$psfile,"-JM -R -Sm0.20 -N", at cmt_files);
+
+plot_pstext(\*BASH,$psfile,"-JM -R -B -G250 -S1p,0","/data2/Datalib/SC/socal_topo_labs.xyz");
+
+plot_psxy(\*BASH,$psfile,"-JX -O","");
+close(BASH);
+system("bash socal_map.bash; gv socal_map.ps");
Property changes on: seismo/3D/GRD_CMT3D/scripts/plot_socal_events.pl
___________________________________________________________________
Name: svn:executable
+ *
Modified: seismo/3D/GRD_CMT3D/scripts/process_data_and_syn.pl
===================================================================
--- seismo/3D/GRD_CMT3D/scripts/process_data_and_syn.pl 2008-12-10 13:50:34 UTC (rev 13642)
+++ seismo/3D/GRD_CMT3D/scripts/process_data_and_syn.pl 2008-12-10 21:07:27 UTC (rev 13643)
@@ -1,5 +1,6 @@
#!/usr/bin/perl -w
+# again this is partly based on Carl's master script
use Getopt::Std;
use POSIX;
use vars qw($opt_0 $opt_1 $opt_2 $opt_3 $opt_t $opt_x);
More information about the CIG-COMMITS
mailing list