velestuser.ektext

这是一个不错的源码

the total possible readings. Knowledge about the near-surface velocity 

structure beneath the reference station may help to qualitatively inter-

prete the station delays.



ATTENTION: For reasons of setup of the ray tracer, all stations MUST lie 

within layer 1 and MUST NOT be at same depth as top of layer 2. Thus, if you 

do not use station elevations put top of layer 2 to +0.05 km in the model file.





6.3  Content of velocity model file ( *.MOD)



The model file consists of a title line and for each model of one line with 

the number of layers (layer boundaries) and for each layer one line with 

the velocity (in km/s), the upper interface of the velocity layer (in km), 

vdamp, reflchr, titl. REFLCHR denotes the reflected phase type, currently 

only used for single_event_mode.



Example of model file:



GABILANS1D-model (mod.06 CT070194)   Ref. station BCGM      

 10        vel,depth,vdamp     (f5.2,5x,f7.2,2x,f7.3,3x,a1)

 4.80	-3.0	001.00	P-VELOCITY MODEL

 5.00	0.05	001.00

 5.20	1.0	001.00

 5.40	3.0	001.00

 5.60	5.0	001.00

 5.90	7.0	001.00

 6.10	10.0	001.00

 6.50	15.0	001.00

 6.80	20.0	001.00

 8.05	28.0	001.00

 10        vel,depth,vdamp     (f5.2,5x,f7.2,2x,f7.3,3x,a1)

 2.80	-3.0	001.00	S-VELOCITY MODEL

 3.00	0.05	001.00

 3.20	1.0	001.00

 3.40	3.0	001.00

 3.60	5.0	001.00

 3.90	7.0	001.00

 4.10	10.0	001.00

 4.50	15.0	001.00

 4.80	20.0	001.00

 4.05	28.0	001.00



Separate damping factors (vdamp) may be applied to each velocity layer. 

vdamp is a real number between 0 and 999.99.  The effective damping 

factor for each velocity layer is calculated as follows: dampeff = damping 

(VTHET, given in file VELEST.CMN) * vdamp.

Thus, if vdamp = 1 --> Ordinary damping as specified in VELEST.CMN, 

vdamp = 0 -->  No damping is applied to this velocity layer, vdamp = 999 -

-> The velocity layer is almost fixed to the starting value (overdamped).







6.4  Content of local earthquake file ( *.CNV)



The input-data in *.CNV - format contains a summary-card of the format: 



(3i2.2,1x,2i2.2,1x,f5.2,1x,f7.4,a1,1x,f8.4,a1,f7.2,f7.2,i2...)

iyear,month,iday,ihr,min,sec,xlat,cns,xlon,cew,depth,emag,ifx...

	(where cns is 'N' or 'S' and cew is 'E' or 'W')



This summary card is followed by an unspecified number of lines each 

with six phase data. The end of the event is marked by a blank line. For 

each phase the station name, the phase type, the observation weight, and 

the observation time (in sec) relative to the origin time must be specified.



Example of earthquake data file in *.CNV format:



90 1 1 1349  0.36 36.6745N 121.3038W    6.08   2.10       1                                     

BLRMP0	1.51	BCGMP0	1.59	BSCMP0	1.62	BHRMP0	2.42	BJOMP0	1.71	BEHMP0	3.22

BSLMP0	3.33	BVYMP0	2.59	BVLMP0	3.08	HJSMP1	3.88 	BJCMP0	2.91	BEMMP2	4.19

HFPMP2	3.32	BSRMP0	3.43	HQRMP1	4.55	HLTMP1	5.39	BPIMP0	4.27	BAVMP0	5.36



90 1 2  310  7.67 36.6922N 121.3207W    4.32   3.20        0                                     

BCGMP0	1.11	BLRMP0	1.58	BHRMP0	2.23	BSCMP0	1.85	BJOMP0	2.37	BSLMP0	2.88





Switch  IFX:   Default value=0



For each event the switch IFX (i2, after magn.) may be set. If IFX=1, then 

adjustement of this hypocenter in y-direction (to rotate y-axis see parameter 

ROTATE in control file *.CMN) is inhibited.



In the same event only one phase of the same type (P or S) is accepted for 

the same station!

VELEST may also use S-P travel time differences. They must be written 

into a phase file of INPUTSED format (ISED=2). The phase description is "S-

P". The phase identifier is a "-" rather than P or S.

Reflected waves are computed only if the arrival-time is specified as one 

of a reflected wave (see below); in this case, the ray is forced to be a 

reflected one and no other raypaths are computed (because a reflected 

wave hardly is the first arrival...!).





7.  Overview over main print OUTPUT file



The head of this file recapitulates the content of the control input file (*.cmn) for 

this run and summarizes the most important parameteres of the short distance 

conversion. Subsequently, the content of the input model file, the station file, and 

the input phase data are listed. This is followed by the so-called 'iteration 0' - a 

one-step forward (location only) solution for each event to obtain the initial 

observation residuals and RMS residuals.

Each iteration starts with a summary of the control parameters and the effective 

number of unknowns. 



number of unknowns = 4*ieq + #station + #blasts + #layers



The number of unknowns might vary between iterations if INVERTRATIO is 

any number other than 1. Before the resulting RMS values and other statistical 

estimators are given a check is performed if the data variance has decreased. If 

new iteration leads to a worse solution than the previous one, model and 

hypocenter are 'backuped'. 

For each event the hypocenter adjustements are listed and at the end the 

average and the average absolute adjustement are calculated to allow an 

overview over the the step length.

The average velocity calculations may provide a summary of the model for 

comparison with known average crustal velocities. Since the first layer 

contains all stations (at possibly variable elevation) and in some cases this 

layer may be very thin and of very low velocity, a second list with average 

velocities to different depth is provided beginning with only layer two.

The information for each iteration ends with the list of the new station 

delays and the adjustements to the old station delays.



References:



Crosson, R.S.: Crustal structure modeling of earthquake data,1, Simultaneous 

least squares estimation of hypocenter and velocity parameters. J. 

Geophys. Res., 81, 3036-3046,1976.

Ellsworth, W.L.: Three-dimensional structure of the crust and mantle beneath the 

island of Hawaii. Ph D thesis, MIT, Massachussetts, USA,1977.

Kissling, E.: Geotomography with local earthquake data, Rev. of Geophysics, 

1988.

Kissling, E., Ellsworth, W.L., and R. Cockerham: Three-dimensional structure 

of the Long Valley Caldera, California, region by geotomography. 

U.S. Geol. Surv. Open File Rep. 84-939, 188-220,1984.

Kissling, E., and J. Lahr: Tomographic image of the Pacific slab under southern 

Alaska. Eclogae Geol. Helv., 84, 297-315, 1991.

Kissling, E., Ellsworth, W. L., Eberhart-Phillips, D., and U. Kradolfer: Initial 

reference models in local earthquake tomography, J. Geophys. Res., 

99, 19'635-19'646, 1994.

Kradolfer, U.: Seismische Tomographie in der Schweiz mittels lokaler Erdbeben. 

Ph D thesis, ETH Zuerich, Switzerland, 1989.

Maurer, H.: Seismotectonics and upper crustal structure in the western Swiss 

Alps. Ph D thesis, ETH Zuerich, Switzerland, 1993.

Roecker, S.W.: Seismicity and tectonics of the Pamir-Hindu Kush region of 

central Asia. Ph D thesis, MIT, Massachussetts, USA,1977.

Thurber, C.H.: Earth structure and earhtquake locations in the Coyote Lake area, 

central California. Ph D thesis, MIT, Massachussetts, USA,1977.







8.  Explanation of SINGLE-EVENT-MODE





Modified to use first station (+.001 deg lat/lon) as trial hypococenter 

(optionally) instead of the one on the first line of each event in the *.CNV 

input-file. To do so, set ITRIAL to 1 and ZTRIAL to the desired value, else 

set ITRIAL to zero.

In case of the hypocenter-parameters, the depth must not be reset by 1/2 

of the step-vector ( b(j) ), but rather by 1/2 of the effective z-change. In 

case of a 'depth-constraint',the depth used to be put down to much during 

the backup-procedure. Now the effective depth-change is stored in 

variable EFFDELTAZ(ievent).

Instead of calculating a EFFDELTAZ in subr. adjustmodel, the depth is in all 

cases constrained directly on the adjustment-vector-element. So in subr. 

backup the hypocentral depth is 'backed up' in the same way as are all the 

other hypocentral parameters. As a consequence of doing so, the calculated 

step-length reflects now the ACTUAL (applied) step-length. Therefore 

further iterations are avoided, if the epicenter is already o.k. but the 

inversion-process continues trying to put the hypocenter above ZMIN . 

Now the step-length gets indeed small if no more adjustments are allowed 

and the iteration process stops as desired in such a case! 



 outputfile (unit=2) with STATISL-compatible output is generated 

(STATISL is a statistic-analysis-program for locations of local earthquake-

data). 

Resolution- and covariance-matrices are calculated and printed in output-

file. Dirichlet-spread of resolution-matrixis determined and added to the 

output (file 16 & screen). 

During first 2 iterations, depth-adjustment is set to 0.0. The same is done if 

igap of the event is > 250 degrees. If depth-adjustment > zadj (set in input-

file) then depth-adjustment is constrained ! 

Covariance matrix is now calculated as the UNIT covariance matrix, so the 

resulting SIGMAs (for OT, X, Y, Z) need to be multiplied by the effective 

(assumed) SQRT(data variance) SIGMAdata (e.g. 0.1) ! Before, cov was 

multiplied by davar1, which occured to be very small sometimes (in case, 

the model-parameters fitted the model well) and therefore unrealistic 

small values for COV were obtained. 

For Swiss-data the local magnitude is calculated (see subr. MAGNITUDE 

and others). Ml is calculated only in single-event-mode 

In single-event-mode the regionname is determined according the 

coordinates of the final solution. Usually the Flinn-Engdahl names are used; 

if possible, the names of the 1:25'000 maps of Switzerland are used. 



In single-event-mode only readings with obs-weight < 4 are used for the 

solution and only those are used for AVRES, RMS, DATVAR and so on. 

Therefore the normalization of the weights is not based on knobs directly, 

but rather on (knobs-nobswithw0), where nobswithw0 is the number of 

observations with w(nobs,ievent)=0.0 [these observations have either an 

observation-weight of 4 or they were rejected during the iteration process 

(see subr. rejectobs). However, these readings are taken for calculation of 

residuals and magnitude. 

File 2 contains first a STATISL-compatible output and second this output is 

equivalent to the one from HYPO71 (Version at the Swiss Seismological 

Service, ETH-Zurich). Two additional parameters are printed out: The ALE 

(Average logarithmic eigenvalue of GtG) and the D-SPR (Dirichlet

cu spread). 

Reject observation if NITT > 2 .and. RES(nobs) > 2.0 sec.  If a rejected 

arrival gets a residual smaller than 1 sec, then its weight is 'revived' and 

again used in further calculation. 

If any of the hypocentral parameters is fixed, both the ALE and the 

Dirichlet-spread are calculated for the SOLVED parameters only! In file02 

[HYPO71-compatible and STATISL-compatible output] a notification made 

if any of the hypocentral parameter is fixed: the word DEPTH or the words 

LAT, LON and DETH respectively, appear as *DEPTH*, *LAT* and *LON* if 

they are fixed.

One single card of a phase-list consists usually not only of the P-phase, but 

also of the S-phase and sometimes amplitude&period-data. In the case, 

that only an S-arrival is read from the seismogram, many location- 

programs need a 'dummy-P-arrival-time' and the operator then usually 

sets the fantasy-P-arrival to an observation-weight of 4! Program VELEST 

allows, to add only a S-arrival (plus eventually amplitude- and period-

data) and NO P-arrival-time. In this case the field, where usually the P-

observation-weight stands, the number  5  must be typed in the phaselist. 

If done so, all the data on the phaselist-line are used with the exception of 

the P-time. Of course, in the output-file (file02) there isn't printed any 

information about any P of this station, but the full information of the S-

ray is printed out. If the input-line contains only P-data and its 

observation-weight is set to 5, the observation isn't used at all.   

  The user of VELEST is recommended (in case the subroutine which inputs 

the phase-data - here INPUTSED - is substitutet by another one) to be sure 

that in case the field with the P-arrival-second is empty, this 'phase' is 

either skipped or set to a weight of 5 !

 If a trial-hypocenter for event-location is to be found, the first station in 

the data which is on the station-list and has an observation-weight less 

than 5 is taken as starting point. Now all subsequent arrivals

 are checked whether they arrive earlier and have an obs.-weight less than 

5 and are on the station-list. The station with the first arrival in time 

which satisfies these conditions is taken as trial epicenter. Note: In case 

that anything goes wrong during this search and the first station is 

recognized to be blank, the first station in the station-list is taken as trial 

epicenter (well, this should almost never be the case, but one never 

knows...).  



Single-event-mode: If depth is fixed and fix-depth is less/equal zero, 

depth is always set to ZMIN (=zmin allowed). ZMIN is read from the input- 

command-file and altered (during iteration-process) to the real surface-

height (depending on current horizointal coordinates) if option ITOPO is 

set. 

The calculation of the data resolution matrix in the single event mode was 

applied. The diagonal elements (importance) are applied to the output of 

the statisl file. The importance is a measurement of the 'importance' of a 

station residual. Large residuals with a high importance are very 

suspicious!!! (See PhD thesis H. Maurer 1993, for further information.) The 

dataresolution (hat) matrix is set up as follows: The matrix G (data kernel) 

is recreated in vel_main.f. Then the damped normal equations are set up. 

Next, the normal equations are solved using the matrix inversion routine 

matrinv (declared in matrinv.f). Finally the hat matrix is calculated and the 

diagonal elements are extracted. The importance values are written into 

the statisl output file (last parameter of each 'phase line').