syndat.f -- sum greens functions with moment tensor and centroid time correction
click title to see Google Drive directory JParkCodes
This code is an updated version of a code first posted on 10/30/2009. The differences are mainly cosmetic, but have practical value. The revision allows the code to circumvent an annoying convention in horizontal components and spherical-earth mode-sum synthetic seismograms. Because the code grn_sph.f writes greensfunctions for radial and transverse horizontal components of motion, it needs to apply the instrument response for one of the actual components of the station. The code was originally written to use the E-component response file. A snag arises when there is no east component for the station because the horizontals are deflected too far from the compass orientations. If the RESP-file for the E-component does not exist, syndat looks for the response file of the "1" component instead. (A station may collect data on the LH1,LH2 components, instead of the LHE,LHN components.) This workaround seems to work in test examples.
ALSO, the code is using the software library from evalresp3.3.3, obtained from IRIS. The C-routine evresp.c in that package spits out a comment that "velocity (m/s)' is being converted to instrument counts. This comment is not correct -- the syndat code tells the evresp() function that its input is displacement, and the application of instrument response looks right. (I also tested the false specification of "VEL" input, but that choice misfit seismic data greatly.
This code is a version of the old sac_syndat.f that incorporates the IRIS-suppied EVRESP subroutine to compute and apply instrument response correction to convert displacement to instrument counts.
c f77 -o /Users/jjpark/bin/syndat syndat.f /Users/jjpark/Ritz/libevresp.a /Users/jjpark/Ritz/jlib.a
c xf77 -o /Users/jjpark/bin/syndat syndat.f /Users/jjpark/Ritz/libevresp.a /Users/jjpark/Plotxy/plotlib.a /Users/jjpark/Ritz/jlib.a
c program to make synthetics using green fcts
c the greenfcts are raw displacement in meters. This code applies the instrument
c response in the frequency domain with a FFT --> multiply by H(f) --> IFFT process.
program asks for the centroid time delay, which is applied with a phase-ramp in the frequency domain
The centroid time delay occurs because the start time of a large event is typically not the best time to place the stepfunction onset of the long-period course. If moment release occurs over T seconds, then a better fit to data is found by delaying the synthetic onset by T/2 seconds. (Assuming a symmetric source time function, natch).
print *,'enter time delay e.g. centroid-time correction'
read(5,*) tconst
40 print *,'do you want to use strike,dip and rake (0=no) '
read(5,*) ick
if(ick.ne.0) then
call fault(sol)
else
print *,'enter f(x6) : '
read(5,*) (sol(i),i=1,6)
endif
print *,'routine expects the gfct filename to start with g'
print *,'the code substitutes s for g in the output filename'
30 print *,'green fns file name(=stop) : '
The program will read a succession of greensfct files, one per station-component, until it reads 'stop' e.g.
0 --> zero centroid time delay
1 --> read strike,dip,slip
329 8 110 --> strike dip slip
400 --> seismic moment (units of 10**27 dyne-cm/10**20 nt-m)
/Users/jjpark/Synth/gsph_ffc.lhz --> greensfct file
/Users/jjpark/Synth/gsph_ffc.lhr
/Users/jjpark/Synth/gsph_ffc.lht
stop
There are PLOTIT calls within the code (commented out for now) that allow the user to view his/her handiwork, testing the causality of the instrument response, etc.
The EVRESP routine is found in the libevresp.a library that you can create from the IRIS EVALRESP software distribution. For reasons only the DMC programmers know, the latest versions of the EVALRESP distribution hide the libevresp.a library in a subdirectory ".lib" that wont show up in a "ls" command unless you are looking for it. Just another convenience feature, I guess.
This code is an updated version of a code first posted on 10/30/2009. The differences are mainly cosmetic, but have practical value. The revision allows the code to circumvent an annoying convention in horizontal components and spherical-earth mode-sum synthetic seismograms. Because the code grn_sph.f writes greensfunctions for radial and transverse horizontal components of motion, it needs to apply the instrument response for one of the actual components of the station. The code was originally written to use the E-component response file. A snag arises when there is no east component for the station because the horizontals are deflected too far from the compass orientations. If the RESP-file for the E-component does not exist, syndat looks for the response file of the "1" component instead. (A station may collect data on the LH1,LH2 components, instead of the LHE,LHN components.) This workaround seems to work in test examples.
ALSO, the code is using the software library from evalresp3.3.3, obtained from IRIS. The C-routine evresp.c in that package spits out a comment that "velocity (m/s)' is being converted to instrument counts. This comment is not correct -- the syndat code tells the evresp() function that its input is displacement, and the application of instrument response looks right. (I also tested the false specification of "VEL" input, but that choice misfit seismic data greatly.
This code is a version of the old sac_syndat.f that incorporates the IRIS-suppied EVRESP subroutine to compute and apply instrument response correction to convert displacement to instrument counts.
c f77 -o /Users/jjpark/bin/syndat syndat.f /Users/jjpark/Ritz/libevresp.a /Users/jjpark/Ritz/jlib.a
c xf77 -o /Users/jjpark/bin/syndat syndat.f /Users/jjpark/Ritz/libevresp.a /Users/jjpark/Plotxy/plotlib.a /Users/jjpark/Ritz/jlib.a
c program to make synthetics using green fcts
c the greenfcts are raw displacement in meters. This code applies the instrument
c response in the frequency domain with a FFT --> multiply by H(f) --> IFFT process.
program asks for the centroid time delay, which is applied with a phase-ramp in the frequency domain
The centroid time delay occurs because the start time of a large event is typically not the best time to place the stepfunction onset of the long-period course. If moment release occurs over T seconds, then a better fit to data is found by delaying the synthetic onset by T/2 seconds. (Assuming a symmetric source time function, natch).
print *,'enter time delay e.g. centroid-time correction'
read(5,*) tconst
40 print *,'do you want to use strike,dip and rake (0=no) '
read(5,*) ick
if(ick.ne.0) then
call fault(sol)
else
print *,'enter f(x6) : '
read(5,*) (sol(i),i=1,6)
endif
print *,'routine expects the gfct filename to start with g'
print *,'the code substitutes s for g in the output filename'
30 print *,'green fns file name(
The program will read a succession of greensfct files, one per station-component, until it reads 'stop' e.g.
0 --> zero centroid time delay
1 --> read strike,dip,slip
329 8 110 --> strike dip slip
400 --> seismic moment (units of 10**27 dyne-cm/10**20 nt-m)
/Users/jjpark/Synth/gsph_ffc.lhz --> greensfct file
/Users/jjpark/Synth/gsph_ffc.lhr
/Users/jjpark/Synth/gsph_ffc.lht
stop
There are PLOTIT calls within the code (commented out for now) that allow the user to view his/her handiwork, testing the causality of the instrument response, etc.
The EVRESP routine is found in the libevresp.a library that you can create from the IRIS EVALRESP software distribution. For reasons only the DMC programmers know, the latest versions of the EVALRESP distribution hide the libevresp.a library in a subdirectory ".lib" that wont show up in a "ls" command unless you are looking for it. Just another convenience feature, I guess.
No comments:
Post a Comment