Friday, October 30, 2009

syndat.f -- sum greens functions with moment tensor and centroid time correction

click title to see Google Drive directory JParkCodes

NOTE:  This code was updated on 7/24/2012.  The title of the post points to the newer code.  The older code can be found at Google Drive JParkCodes/syndat_2009.f

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.

No comments:

 
Link