Non-linear hypocenter location with UTM geographical projection. Used for locating earthquakes in local or regional settings.

1 2 3 4 5 |

`Ldat` |
list, must inlude: x,y,err, sec, cor (see details) |

`EQ` |
list, must inlude: x,y,z, t |

`vel` |
list, 1D velocity structure |

`maxITER` |
Maximum number of iterations |

`distwt` |
distance weighting factor |

`lambdareg` |
regularization parameter for damping |

`FIXZ` |
logical, TRUE = fix depth, i.e. only calculate x,y,t |

`REG` |
logical, TRUE=use regularization |

`WTS` |
logical, TRUE==use weighting |

`STOPPING` |
logical, TRUE=use stopping criteria |

`RESMAX` |
vector, residual max for P and S, default=c(4,5) |

`tolx` |
numeric, tolerance in km in x direction |

`toly` |
numeric, tolerance in km in y direction |

`tolz` |
numeric, tolerance in km in z direction |

`PLOT` |
logical, plot results during iterations |

Input pick list must have at x,y,z, sec, cor, err elements for each station. If no station correction is available it is set to zero. If no uncertainty (err) is available, it is set to 0.05 sec. Each station must have a finite x-y coordinate and arrival time in seconds. Events are located relative to the minute.

Routine uses the svd in a sequence of linear inversions to estimate the nonlinear location.

List:

`EQ` |
list, Earthquake hypocenter and time |

`its` |
number of iterations |

`rms` |
rms residual |

`wrms` |
wheighted rms residual |

`used` |
vector, index of used equations |

`guesses` |
list of x,y,z,t intermediate locations when converging |

This routine should be called by a wrapper (Vlocate) that applies the algorithm several times and changes parameters based on the quality.

If RESMAX is used and the robust approach yields fewer than 4 equations, the best (smallest) four residuals will be used to determiine the event location.

Jonathan M. Lees<jonathan.lees@unc.edu>

Vlocate

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 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 | ```
## Not run:
library(RSEIS)
data(GH)
g1 = GH$pickfile
data(VELMOD1D)
vel= VELMOD1D
w1 = which(!is.na(g1$STAS$lat))
sec = g1$STAS$sec[w1]
N = length(sec)
Ldat = list(
name = g1$STAS$name[w1],
sec = g1$STAS$sec[w1],
phase = g1$STAS$phase[w1],
lat=g1$STAS$lat[w1],
lon = g1$STAS$lon[w1],
z = g1$STAS$z[w1],
err= g1$STAS$err[w1],
yr = rep(g1$LOC$yr , times=N),
jd = rep(g1$LOC$jd, times=N),
mo = rep(g1$LOC$mo, times=N),
dom = rep(g1$LOC$dom, times=N),
hr =rep( g1$LOC$hr, times=N),
mi = rep(g1$LOC$mi, times=N) )
MLAT = median(Ldat$lat)
MLON = median(Ldat$lon)
proj = GEOmap::setPROJ(type=2, LAT0=MLAT, LON0=MLON)
#### get station X-Y values in km
XY = GEOmap::GLOB.XY(Ldat$lat, Ldat$lon, proj)
### add to Ldat list
Ldat$x = XY$x
Ldat$y = XY$y
wstart = which.min(Ldat$sec)
EQ = list(x=XY$x[wstart], y=XY$y[wstart], z=6, t=Ldat$sec[wstart] )
maxITER = 7
###print(EQ)
AQ = XYlocate(Ldat,EQ,vel,
maxITER = maxITER,
distwt = 1,
lambdareg =10 ,
FIXZ = FALSE,
REG = TRUE,
WTS = TRUE,
STOPPING = TRUE,
RESMAX = c(0.1,0.1),
tolx = 0.001,
toly = 0.001 ,
tolz = 0.5, PLOT=FALSE)
######## update the new location
AXY = GEOmap::XY.GLOB(AQ$EQ$x, AQ$EQ$y, proj)
AQ$EQ$lat = AXY$lat
AQ$EQ$lon = AXY$lon
if(AQ$EQ$lon>180) { AQ$EQ$lon = AQ$EQ$lon-360 }
plot(c(Ldat$x, AQ$EQ$x) , c(Ldat$y,AQ$EQ$y), type='n' , xlab="km",
ylab="km" )
points(Ldat$x, Ldat$y, pch=6)
points(AQ$EQ$x, AQ$EQ$y, pch=8, col='red')
points(EQ$x, EQ$y, pch=4, col='blue')
legend("topright", pch=c(8,4, 6), col=c("red", "blue", "black"),
legend=c("Final location", "Initial guess", "Station"))
print(AQ)
EQ$x = 10
EQ$y = 2
## End(Not run)
``` |

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.

All documentation is copyright its authors; we didn't write any of that.