[cig-commits] r12791 - seismo/3D/ADJOINT_TOMO/iterate_adj/matlab
carltape at geodynamics.org
carltape at geodynamics.org
Wed Sep 3 11:34:18 PDT 2008
Author: carltape
Date: 2008-09-03 11:34:18 -0700 (Wed, 03 Sep 2008)
New Revision: 12791
Added:
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/OUTPUT_MISFIT/
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/OUTPUT_SUBSPACE/
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/defval.m
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/nounder.m
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/osdep.m
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/readsac.m
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_misfit_dir.pl
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_subspace_dir.pl
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/suf.m
Removed:
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/OUTPUT/
seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/make_dirs.pl
Log:
Changed OUTPUT directory names to prepare for generating comparison plots of data and synthetics for two different models. Also added matlab scripts from Frederik Simons webpage to allow for reading sac files into Matlab -- readsac.m, suf.m, defval.m, osdep.m, nounder.m.
Copied: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/OUTPUT_SUBSPACE (from rev 12790, seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/OUTPUT)
Added: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/defval.m
===================================================================
--- seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/defval.m (rev 0)
+++ seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/defval.m 2008-09-03 18:34:18 UTC (rev 12791)
@@ -0,0 +1,24 @@
+function defval(name,value)
+% DEFVAL(name,value)
+%
+% Checks if a variable 'name' exists in the caller's workspace,
+% if nonexistent or if empty 'name' is set to value 'value' in caller's workspace
+%
+% Won't work for an unassigned structure variable
+%
+% Last modified by fjsimons-at-alum.mit.edu, 09/12/2007
+
+if ~ischar(name),
+ error('The first argument of defval has to be a string (variable name)');
+end
+
+si=1;
+if evalin('caller',[ 'exist(''' name ''')']);,
+ si=evalin('caller',[ 'isempty(' name ')']);
+end
+if si,
+ assignin('caller',name,value);
+ na=dbstack;
+ %disp(['Default value used in ' na(2).name ': ' name '=' num2str(value(1,:))])
+end
+
Deleted: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/make_dirs.pl
===================================================================
--- seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/make_dirs.pl 2008-09-03 15:04:13 UTC (rev 12790)
+++ seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/make_dirs.pl 2008-09-03 18:34:18 UTC (rev 12791)
@@ -1,65 +0,0 @@
-#!/usr/bin/perl -w
-
-#-----------------------------------
-# Carl Tape, 11-July-2008
-# make_dirs.pl
-#
-# This script empties the files within the mu directories.
-# These files are primarily the mu vectors used in the subspace method.
-#
-# This should be run prior to running subspace_specfem.m.
-#
-# EXAMPLE: make_dirs.pl 10
-#
-#-----------------------------------
-
-if (@ARGV < 1) {die("Usage: make_dirs.pl smodel\n")}
-($maxm) = @ARGV;
-
-# possible labels for choice of data covariance
- at dcov_tags = ("event","window","none");
-
-# possible labels for choice of kernels
- at kers = ("mu_kernel","kappa_kernel","mu_kernel_smooth","kappa_kernel_smooth","both","both_smooth");
-
-$ndcov = @dcov_tags;
-$nker = @kers;
-
-$dir0 = "OUTPUT";
-`mkdir -p $dir0`;
-#if(not -e $dir0) {die("dir dir0 $dir0 does not exist\n");}
-
-#$dir1 = "$dir0/$smodel";
-#`mkdir -p $dir1`;
-
-for ($h = 0; $h <= $maxm; $h = $h+1) {
-
- $smodel = sprintf("m%2.2i",$h);
- $dir1 = "${dir0}/${smodel}";
- `mkdir -p $dir1`;
-
- for ($j = 1; $j <= $ndcov; $j = $j+1) {
-
- $dcov = $dcov_tags[$j-1];
- $dir2 = "${dir1}/${dcov}";
- `mkdir -p $dir2`;
-
- for ($i = 1; $i <= $nker; $i = $i+1) {
-
- $ker = $kers[$i-1];
- $dir = "$dir2/${ker}";
- print "$dir \n";
-
- if (-e $dir) {
- #print "--> dir exists -- now deleting and remaking\n"; `rm -rf $dir`;
- print "--> dir exists\n";
- } else {
- print "--> dir does not exist -- now making\n";
- }
- `mkdir -p $dir`;
- `mkdir -p $dir/mu_all`;
- }
- }
-}
-
-#================================================
Added: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/nounder.m
===================================================================
--- seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/nounder.m (rev 0)
+++ seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/nounder.m 2008-09-03 18:34:18 UTC (rev 12791)
@@ -0,0 +1,9 @@
+function sin=nounder(sin)
+% sin=NOUNDER(sin)
+%
+% Changes underscores to dashes in a string
+%
+% Last modified by fjsimons-at-alum.mit.edu, Feb 06th, 2004
+
+sin(find(abs(sin)==95))='-';
+
Added: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/osdep.m
===================================================================
--- seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/osdep.m (rev 0)
+++ seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/osdep.m 2008-09-03 18:34:18 UTC (rev 12791)
@@ -0,0 +1,16 @@
+function of=osdep
+% of=OSDEP
+%
+% Returns the value of the read option to be used
+% on the LOCAL operating system for files created
+% on the LOCAL operating system.
+%
+% Last modified by fjsimons-at-alum.mit.edu, 23.11.2004
+
+
+if strcmp(getenv('OSTYPE'),'linux')
+ of= 'l';
+end
+if strcmp(getenv('OSTYPE'),'solaris')
+ of= 'b';
+end
Added: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/readsac.m
===================================================================
--- seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/readsac.m (rev 0)
+++ seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/readsac.m 2008-09-03 18:34:18 UTC (rev 12791)
@@ -0,0 +1,110 @@
+function varargout=readsac(filename,plotornot,osd)
+% [SeisData,HdrData,tnu,pobj,tims]=READSAC(filename)
+% [SeisData,HdrData,tnu,pobj,tims]=READSAC(filename,1)
+% [SeisData,HdrData,tnu,pobj,tims]=READSAC(filename,1,'l')
+% [SeisData,HdrData,tnu,pobj,tims]=READSAC(filename,1,'b')
+%
+% READSAC(...)
+%
+% Reads SAC-formatted data (and makes a plot)
+%
+% INPUT:
+%
+% filename The filename, full path included
+% plotornot 1 Makes a plot
+% 0 Does not make a plot [default]
+% osd 'b' for data saved on Solaris read into Linux
+% 'l' for data saved on Linux read into Linux
+%
+% OUPUT:
+%
+% SeisData The number vector
+% HdrData The header structure array
+% tnu The handle to the plot title; or appropriate string
+% pobj The handle to the plot line and the xlabel
+% tims The times on the x-axis
+%
+% See SETVAR.PL out of SACLIB.PM
+%
+% Last modified by fjsimons-at-alum.mit.edu, 07/07/2008
+
+defval('plotornot',0)
+defval('osd',osd)
+
+fid=fopen(filename,'r',osd);
+if fid==-1
+ error([ 'File ',filename,' does not exist in current path ',pwd])
+end
+HdrF=fread(fid,70,'float32');
+HdrN=fread(fid,15,'int32');
+HdrI=fread(fid,20,'int32');
+HdrL=fread(fid,5,'int32');
+HdrK=str2mat(fread(fid,[8 24],'char'))';
+SeisData=fread(fid,HdrN(10),'float32');
+fclose(fid);
+
+% Meaning of LCALDA from http://www.llnl.gov/sac/SAC_Manuals/FileFormatPt2.html
+% LCALDA is TRUE if DIST, AZ, BAZ, and GCARC are to be calculated from
+% station and event coordinates. But I don't know how to set this
+% variable to a logical yet using setvar.pl. What is the TRUE/FALSE
+% format in SAC? 18.1.2006
+
+% If you change any of this, change it in WRITESAC as well!
+HdrData=struct(...
+ 'AZ',HdrF(52),...
+ 'B',HdrF(6),...
+ 'BAZ',HdrF(53),...
+ 'DIST',HdrF(51),...
+ 'DELTA',HdrF(1),...
+ 'E',HdrF(7),...
+ 'EVDP',HdrF(39),...
+ 'EVEL',HdrF(38),...
+ 'EVLA',HdrF(36),...
+ 'EVLO',HdrF(37),...
+ 'GCARC',HdrF(54),...
+ 'IFTYPE',HdrI(1),...
+ 'KCMPNM',HdrK(21,HdrK(21,:)>64),...
+ 'KINST',HdrK(24,HdrK(24,:)>64),...
+ 'KSTNM',HdrK(1,HdrK(1,:)>64),...
+ 'KUSER0',HdrK(18,HdrK(18,:)>64),...
+ 'LCALDA',HdrL(4),...
+ 'MAG',HdrF(40),...
+ 'NPTS',HdrN(10),...
+ 'NVHDR',HdrN(7),...
+ 'NZHOUR',HdrN(3),...
+ 'NZJDAY',HdrN(2),...
+ 'NZMIN',HdrN(4),...
+ 'NZMSEC',HdrN(6),...
+ 'NZSEC',HdrN(5),...
+ 'NZYEAR',HdrN(1),...
+ 'SCALE',HdrF(4),...
+ 'STDP',HdrF(35),...
+ 'STEL',HdrF(34),...
+ 'STLA',HdrF(32),...
+ 'STLO',HdrF(33),...
+ 'T0',HdrF(11),...
+ 'T1',HdrF(12),...
+ 'T2',HdrF(13),...
+ 'T3',HdrF(14));
+
+tims=linspace(HdrData.B,HdrData.E,HdrData.NPTS);
+
+if plotornot==1
+ pobj=plot(tims,SeisData);
+ filename(find(abs(filename)==95))='-';
+ tito=nounder(suf(filename,'/'));
+ if isempty(tito)
+ tito=nounder(filename);
+ end
+ tnu=title(tito);
+ axis tight
+ pobj(2)=xlabel([ 'Time (s)']);
+else
+ tnu=nounder(suf(filename,'/'));
+ pobj=0;
+end
+
+nam={'SeisData' 'HdrData' 'tnu' 'pobj','tims'};
+for index=1:nargout
+ varargout{index}=eval(nam{index});
+end
Added: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_misfit_dir.pl
===================================================================
--- seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_misfit_dir.pl (rev 0)
+++ seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_misfit_dir.pl 2008-09-03 18:34:18 UTC (rev 12791)
@@ -0,0 +1,21 @@
+#!/usr/bin/perl -w
+
+#-----------------------------------
+# Carl Tape, 03-Sept-2008
+# setup_misfit_dir.pl
+#
+# This script xxx
+#
+# EXAMPLE: setup_misfit_dir.pl m00 m07 2/40
+#
+#-----------------------------------
+
+if (@ARGV < 3) {die("Usage: setup_misfit_dir.pl smodel_min smodel_max Tmin/Tmax\n")}
+($minm,$maxm,$Tpers) = @ARGV;
+
+$dir0 = "OUTPUT_MISFIT";
+`mkdir -p $dir0`;
+
+# TO DO --->
+
+#================================================
Property changes on: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_misfit_dir.pl
___________________________________________________________________
Name: svn:executable
+ *
Copied: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_subspace_dir.pl (from rev 12790, seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/make_dirs.pl)
===================================================================
--- seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_subspace_dir.pl (rev 0)
+++ seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/setup_subspace_dir.pl 2008-09-03 18:34:18 UTC (rev 12791)
@@ -0,0 +1,70 @@
+#!/usr/bin/perl -w
+
+#-----------------------------------
+# Carl Tape, 03-Sept-2008
+# setup_subspace_dir.pl
+#
+# This script creates a set of directories for computing misfit
+# functions and mu-vectors for the subspace method for different
+# choices of the data covariance matrix.
+#
+# The user should keep a copy of the output files in OUTPUT_SUBSPACE
+# somewhere outside of the working SVN directory, since these output
+# files will not be checked into SVN.
+#
+# This should be run prior to running subspace_specfem.m.
+#
+# EXAMPLE: setup_subspace_dir.pl 10
+#
+#-----------------------------------
+
+if (@ARGV < 1) {die("Usage: setup_subspace_dir.pl smodel\n")}
+($maxm) = @ARGV;
+
+# possible labels for choice of data covariance
+ at dcov_tags = ("event","window","none");
+
+# possible labels for choice of kernels
+ at kers = ("mu_kernel","kappa_kernel","mu_kernel_smooth","kappa_kernel_smooth","both","both_smooth");
+
+$ndcov = @dcov_tags;
+$nker = @kers;
+
+$dir0 = "OUTPUT_SUBSPACE";
+`mkdir -p $dir0`;
+#if(not -e $dir0) {die("dir dir0 $dir0 does not exist\n");}
+
+#$dir1 = "$dir0/$smodel";
+#`mkdir -p $dir1`;
+
+for ($h = 0; $h <= $maxm; $h = $h+1) {
+
+ $smodel = sprintf("m%2.2i",$h);
+ $dir1 = "${dir0}/${smodel}";
+ `mkdir -p $dir1`;
+
+ for ($j = 1; $j <= $ndcov; $j = $j+1) {
+
+ $dcov = $dcov_tags[$j-1];
+ $dir2 = "${dir1}/${dcov}";
+ `mkdir -p $dir2`;
+
+ for ($i = 1; $i <= $nker; $i = $i+1) {
+
+ $ker = $kers[$i-1];
+ $dir = "$dir2/${ker}";
+ print "$dir \n";
+
+ if (-e $dir) {
+ #print "--> dir exists -- now deleting and remaking\n"; `rm -rf $dir`;
+ print "--> dir exists\n";
+ } else {
+ print "--> dir does not exist -- now making\n";
+ }
+ `mkdir -p $dir`;
+ `mkdir -p $dir/mu_all`;
+ }
+ }
+}
+
+#================================================
Added: seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/suf.m
===================================================================
--- seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/suf.m (rev 0)
+++ seismo/3D/ADJOINT_TOMO/iterate_adj/matlab/suf.m 2008-09-03 18:34:18 UTC (rev 12791)
@@ -0,0 +1,39 @@
+function suffix=suf(string,delim,iter)
+% suffix=SUF(string,delim,iter)
+%
+% Finds the suffix or extension in a string, say a filename.
+%
+% INPUT:
+%
+% string The string to be parsed
+% delim The string that separates prefix from suffix
+% iter 0 Don't iterate
+% 1 Iterate on the suffix to get the last bit
+%
+% See also PREF
+%
+% Last modified by fjsimons-at-alum.mit.edu, 05/06/2008
+
+defval('delim','.');
+defval('iter',1)
+
+ldel=length(delim);
+[prefix,suffix]=strtok(string,delim);
+
+% This added 17.11.2004
+if iscell(suffix)
+ suffix=cell2mat(suffix);
+end
+
+suffix=suffix(ldel+1:end);
+
+if iter==1
+ % This suffix might not be the last one
+ while findstr(suffix,delim)
+ suffix=suf(suffix,delim);
+ end
+end
+
+
+
+
More information about the cig-commits
mailing list