1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28 package edu.sc.seis.TauP;
29
30 import java.io.File;
31 import java.io.FileNotFoundException;
32 import java.io.IOException;
33 import java.io.OptionalDataException;
34 import java.io.StreamCorruptedException;
35 import java.util.Vector;
36 import edu.sc.seis.seisFile.sac.SacTimeSeries;
37
38 /***
39 * Calculate times for phases and set sac headers based on gcarc or dist or
40 * station lat and lon and event lat and lon.
41 *
42 * Note that triplicated phases will cause problems, as there is only one
43 * spot to put a time. An improved method would allow a phase to have several
44 * t#'s associated with it, so that all arrivals could be marked. Currently
45 * however, only the first arrival for a phase name is used.
46 *
47 * Warning: I assume the evdp header has depth in meters unless the -evdpkm
48 * flag is set, in which case I assume kilometers. This may be a problem for
49 * users that improperly use kilometers for the depth units. Due to much
50 * abuse of the sac depth header units I output a warning message if the depth
51 * appears to be in kilometers, ie it is < 1000. This can be safely ignored
52 * if the event really is less than 1000 meters deep.
53 *
54 * @version 1.1.3 Wed Jul 18 15:00:35 GMT 2001
55
56
57
58 * @author H. Philip Crotwell
59 *
60 */
61 public class TauP_SetSac extends TauP_Time {
62
63 protected Vector sacFileNames = new Vector(10);
64
65 protected boolean evdpkm = false;
66
67 public boolean getEvdpkm() {
68 return evdpkm;
69 }
70
71 public void setEvdpkm(boolean evdpkm) {
72 this.evdpkm = evdpkm;
73 }
74
75 public void setSacFileNames(String[] sacFileNames) {
76 this.sacFileNames.removeAllElements();
77 for (int i=0; i<sacFileNames.length; i++) {
78 this.sacFileNames.addElement(sacFileNames[i]);
79 }
80 }
81
82 protected TauP_SetSac() {
83 super();
84 }
85
86 public TauP_SetSac(TauModel tMod) throws TauModelException {
87 super(tMod);
88 }
89
90 public TauP_SetSac(String modelName) throws TauModelException {
91 super(modelName);
92 }
93
94 protected void setSacVarNums() {
95 boolean[] headersUsed = new boolean[10];
96 for (int i=0;i<headersUsed.length;i++) {
97 headersUsed[i] = false;
98 }
99
100 for (int i=0; i<phaseNames.size() && i<10; i++) {
101 if (((PhaseName)phaseNames.elementAt(i)).sacTNum != -1) {
102 headersUsed[((PhaseName)phaseNames.elementAt(i)).sacTNum] = true;
103 }
104 }
105 int j;
106 for (int i=0 ; i<phaseNames.size() ; i++) {
107 if (((PhaseName)phaseNames.elementAt(i)).sacTNum == -1) {
108
109 for (j=0;j<headersUsed.length && headersUsed[j];j++) {}
110 if (j<10) {
111 ((PhaseName)phaseNames.elementAt(i)).sacTNum = j;
112 headersUsed[j] = true;
113 }
114 }
115 }
116
117 }
118
119 public void calculate(double degrees) {
120 recalcPhases();
121 calcTime(degrees);
122 }
123
124 public void init() throws IOException {
125 super.init();
126 setSacVarNums();
127 }
128
129 public void start() throws IOException, TauModelException {
130 SacTimeSeries sacFile = new SacTimeSeries();
131 int phaseNum;
132 double deg;
133
134 for (int i=0;i<sacFileNames.size();i++) {
135 System.out.println((String)sacFileNames.elementAt(i));
136 sacFile.read((String)sacFileNames.elementAt(i));
137 if (sacFile.evdp == -12345.0f) {
138 System.out.println("Depth not set in "+
139 (String)sacFileNames.elementAt(i)+", skipping");
140 continue;
141 }
142 if (sacFile.o == -12345.0f) {
143 System.out.println("O marker not set in "+
144 (String)sacFileNames.elementAt(i)+", skipping");
145 continue;
146 }
147
148 if (sacFile.gcarc != -12345.0f) {
149 if (verbose) {
150 System.out.println("Using gcarc: "+sacFile.gcarc);
151 }
152 deg = sacFile.gcarc;
153 } else if (sacFile.dist != -12345.0f) {
154 if (verbose) {
155 System.out.println("Using dist: "+sacFile.dist);
156 }
157 deg = sacFile.dist/6371.0*180.0/Math.PI;
158 } else if (sacFile.stla != -12345.0f && sacFile.stlo != -12345.0f &&
159 sacFile.evla != -12345.0f && sacFile.evlo != -12345.0f) {
160 if (verbose) {
161 System.out.println("Using stla,stlo, evla,evlo to calculate");
162 }
163 Alert.warning("Warning: Sac header gcarc is not set,",
164 "using lat and lons to calculate distance.");
165 Alert.warning("No ellipticity correction will be applied.",
166 "This may introduce errors. Please see the manual.");
167 deg = SphericalCoords.distance(sacFile.stla, sacFile.stlo,
168 sacFile.evla, sacFile.evlo);
169 } else {
170
171 Alert.warning("Can't get a distance, all distance fields are undef.",
172 "skipping "+(String)sacFileNames.elementAt(i));
173 continue;
174 }
175
176 if (!((evdpkm && depth == sacFile.evdp) ||
177 (!evdpkm && depth == 1000*sacFile.evdp))) {
178 if (!evdpkm && sacFile.evdp != 0 && sacFile.evdp < 1000.0) {
179 Alert.warning("Sac header evdp is < 1000 in "+
180 (String)sacFileNames.elementAt(i),
181 "If the depth is in kilometers instead of meters "+
182 "(default), you should use the -evdpkm flag");
183 }
184 if (evdpkm) {
185 depthCorrect(sacFile.evdp);
186 } else {
187 depthCorrect(sacFile.evdp/1000.0);
188 }
189 }
190 if (verbose) {
191 System.out.println(sacFileNames.elementAt(i)+
192 " searching for "+getPhaseNameString());
193 }
194 calculate(deg);
195
196 if (verbose) {
197 System.out.println(sacFileNames.elementAt(i)+
198 " "+arrivals.size()+" arrivals found.");
199 }
200 for (int arrivalNum=arrivals.size()-1;arrivalNum>=0;arrivalNum--) {
201 phaseNum = -1;
202 for (int j=phaseNames.size()-1;j>=0;j--) {
203 if (getArrival(arrivalNum).name.equals(
204 ((PhaseName)phaseNames.elementAt(j)).name)) {
205 phaseNum = j;
206 break;
207 }
208 }
209 if (phaseNum != -1) {
210 if (verbose) {
211 System.out.println(sacFileNames.elementAt(i)+" phase found "+
212 getArrival(arrivalNum).name+" -> t"+
213 ((PhaseName)phaseNames.elementAt(phaseNum)).sacTNum +
214 ", travel time=" + (float)getArrival(arrivalNum).time);
215 }
216 setSacTHeader(sacFile, ((PhaseName)phaseNames.elementAt(phaseNum)).sacTNum, getArrival(arrivalNum));
217 }
218 }
219 sacFile.write((String)sacFileNames.elementAt(i));
220 }
221 }
222
223 public static void setSacTHeader(SacTimeSeries sacFile, int headerNum, Arrival arrival) {
224 switch(headerNum) {
225 case 0:
226 sacFile.t0 = sacFile.o + (float)arrival.time;
227 sacFile.kt0 = arrival.name;
228 sacFile.user0 = (float)arrival.getRayParam();
229 break;
230 case 1:
231 sacFile.t1 = sacFile.o + (float)arrival.time;
232 sacFile.kt1 = arrival.name;
233 sacFile.user1 = (float)arrival.getRayParam();
234 break;
235 case 2:
236 sacFile.t2 = sacFile.o + (float)arrival.time;
237 sacFile.kt2 = arrival.name;
238 sacFile.user2 = (float)arrival.getRayParam();
239 break;
240 case 3:
241 sacFile.t3 = sacFile.o + (float)arrival.time;
242 sacFile.kt3 = arrival.name;
243 sacFile.user3 = (float)arrival.getRayParam();
244 break;
245 case 4:
246 sacFile.t4 = sacFile.o + (float)arrival.time;
247 sacFile.kt4 = arrival.name;
248 sacFile.user4 = (float)arrival.getRayParam();
249 break;
250 case 5:
251 sacFile.t5 = sacFile.o + (float)arrival.time;
252 sacFile.kt5 = arrival.name;
253 sacFile.user5 = (float)arrival.getRayParam();
254 break;
255 case 6:
256 sacFile.t6 = sacFile.o + (float)arrival.time;
257 sacFile.kt6 = arrival.name;
258 sacFile.user6 = (float)arrival.getRayParam();
259 break;
260 case 7:
261 sacFile.t7 = sacFile.o + (float)arrival.time;
262 sacFile.kt7 = arrival.name;
263 sacFile.user7 = (float)arrival.getRayParam();
264 break;
265 case 8:
266 sacFile.t8 = sacFile.o + (float)arrival.time;
267 sacFile.kt8 = arrival.name;
268 sacFile.user8 = (float)arrival.getRayParam();
269 break;
270 case 9:
271 sacFile.t9 = sacFile.o + (float)arrival.time;
272 sacFile.kt9 = arrival.name;
273 sacFile.user9 = (float)arrival.getRayParam();
274 break;
275 default:
276 break;
277 }
278 }
279
280 public void printStdUsage() {
281 String className = this.getClass().getName();
282 className =
283 className.substring(className.lastIndexOf('.')+1,className.length());
284
285 System.out.println("Usage: "+className.toLowerCase()+" [arguments]");
286 System.out.println(" or, for purists, java "+this.getClass().getName()+
287 " [arguments]");
288 System.out.println("\nArguments are:");
289
290 System.out.println(
291 "-ph phase list -- comma separated phase list,\n"+
292 " use phase-# to specify the sac header,\n"+
293 " for example, ScS-8 puts ScS in t8\n"+
294 "-pf phasefile -- file containing phases\n\n"+
295 "-mod[el] modelname -- use velocity model \"modelname\" for calculations\n"+
296 " Default is iasp91.\n\n");
297 }
298
299 public void printStdUsageTail() {
300 System.out.println(
301 "\n-debug -- enable debugging output\n"+
302 "-verbose -- enable verbose output\n"+
303 "-version -- print the version\n"+
304 "-help -- print this out, but you already know that!\n");
305 }
306
307 public void printUsage() {
308 printStdUsage();
309 System.out.println(
310 "-evdpkm -- sac depth header is in km, default is meters\n");
311 printStdUsageTail();
312 System.out.println(
313 "sacfilename [sacfilename ...]");
314 System.out.println("\nEx: taup_setsac "+
315 "-mod S_prem -ph S-8,ScS-9 wmq.r wmq.t wmq.z");
316 System.out.println("puts the first S arrival in T8 and ScS in T9");
317 }
318
319 public String[] parseCmdLineArgs(String[] args) throws IOException {
320 int i=0;
321 String[] leftOverArgs;
322 int numNoComprendoArgs = 0;
323 File tempFile;
324
325 leftOverArgs = super.parseCmdLineArgs(args);
326 String[] noComprendoArgs = new String[leftOverArgs.length];
327
328 while (i<leftOverArgs.length) {
329 if (leftOverArgs[i].equalsIgnoreCase("-evdpkm")) {
330 evdpkm = true;
331 } else if (leftOverArgs[i].equals("-help")) {
332 noComprendoArgs[numNoComprendoArgs++] = leftOverArgs[i];
333 } else {
334 tempFile = new File(leftOverArgs[i]);
335 if (tempFile.exists() && tempFile.isFile() &&tempFile.canRead()) {
336 sacFileNames.addElement(leftOverArgs[i]);
337 } else {
338 noComprendoArgs[numNoComprendoArgs++] = leftOverArgs[i];
339 }
340 }
341 i++;
342 }
343
344 if (numNoComprendoArgs > 0) {
345 String[] temp = new String[numNoComprendoArgs];
346 System.arraycopy(noComprendoArgs,0,temp,0,numNoComprendoArgs);
347 return temp;
348 } else {
349 return new String[0];
350 }
351 }
352
353 /*** Allows TauP_SetSac to run as an application. Creates an instance
354 * of TauP_SetSac.
355 * . */
356 public static void main(String[] args)
357 throws FileNotFoundException,
358 IOException,
359 StreamCorruptedException,
360 ClassNotFoundException,
361 OptionalDataException
362 {
363 TauP_SetSac tauPSetSac = new TauP_SetSac();
364 if (args.length==0) {
365 tauPSetSac.printUsage();
366 System.exit(1);
367 } else try {
368
369 String[] noComprendoArgs = tauPSetSac.parseCmdLineArgs(args);
370 if (noComprendoArgs.length > 0) {
371 for (int i=0;i<noComprendoArgs.length;i++) {
372 if (noComprendoArgs[i].equals("-help") ||
373 noComprendoArgs[i].equals("-version")) {
374 System.exit(0);
375 }
376 }
377 System.out.println("I don't understand the following arguments, continuing:");
378 for (int i=0;i<noComprendoArgs.length;i++) {
379 System.out.print(noComprendoArgs[i]+" ");
380 if (noComprendoArgs[i].equals("-help") ||
381 noComprendoArgs[i].equals("-version")) {
382 System.out.println();
383 System.exit(0);
384 }
385 }
386 System.out.println();
387 noComprendoArgs = null;
388 }
389
390 if (tauPSetSac.DEBUG) {
391 System.out.println("Done reading "+tauPSetSac.modelName);
392 }
393
394 tauPSetSac.init();
395
396 tauPSetSac.start();
397
398 } catch (TauModelException e) {
399 System.out.println("Caught TauModelException: "+e.getMessage());
400 e.printStackTrace();
401 }
402
403 }
404 }