Rquake (version 2.5-1)

EQXYresid: Calculate Residuals

Description

Given an earthquake hypocenter and a list of station information, retrieve the station residuals.

Usage

EQXYresid(XY, vel = list(), h1 = c(0, 0, 0, 0), PLOT = FALSE)

Value

vector, right hand side of the least squares problem.

Arguments

XY

matrix of station location and arrival times.

vel

list, RSEIS velocity model

h1

hypocenter location, c(x,y,z,t)

PLOT

logical, TRUE=plot the residuals

Author

Jonathan M. Lees<[e-mail]>

Details

The XY mtrix is in cartesian coordinates, i.e. it has been projected into units of km. Only 1D velocity models are used at this time. Only residuals of P and S wave arrivals are estimated.

See Also

travel.time1D,UPdateEQLOC

Examples

####  get sample data
data(GH, package='RSEIS')

pstas  = GH$pickfile

######  get velocity file
v = GH$velfile

####  project to flatten
proj = GEOmap::setPROJ(type = 2, LAT0 = mean(pstas$STAS$lat), LON0 = mean(pstas$STAS$lon) )

 XY = GEOmap::GLOB.XY(pstas$STAS$lat, pstas$STAS$lon, proj)
#######  elevation corrections
    elcor = rep(0, length(pstas$STAS$lat))
    DZ = pstas$STAS$z - mean(pstas$STAS$z)
    elcor[pstas$STAS$phase=="P"] = DZ[pstas$STAS$phase=="P"]/v$vp[1]
    elcor[pstas$STAS$phase=="S"] = DZ[pstas$STAS$phase=="S"]/v$vs[1]

######   set up requisite vectors
    XY$cor = elcor
    XY$phase = pstas$STAS$phase
    XY$sec = pstas$STAS$sec

    sol = c(GH$pickfile$LOC$lat, GH$pickfile$LOC$lon, GH$pickfile$LOC$z, GH$pickfile$LOC$sec)
    
    eqXY = GEOmap::GLOB.XY(sol[1], sol[2], proj)

#######  get residuals
    res =  EQXYresid(XY, vel=v , h1=c(eqXY$x, eqXY$y, sol[3], sol[4] ) ,
    PLOT=FALSE)