Download or view forensicSeismology.frink in plain text format
// Solution to the "Forensic Seismology" geocache:
//
// http://www.geocaching.com/seek/cache_details.aspx?guid=68083e90-f516-46de-a2c9-55709818a7ed
// Import the sun.frink file which contains great circle calculations
// and definitions for many of the functions below.
use navigation.frink
// Coordinates of seismic stations
lata = DMS[39,30.840] North
longa = DMS[104,44.109] West
latb = DMS[39, 31.710] North
longb = DMS[104, 47.901] West
latc = DMS[39, 27.302] North
longc = DMS[104, 46.119] West
// Find the distances between seismic stations
[dab, bab] = earthDistanceAndBearing[lata, longa, latb, longb]
[dbc, bbc] = earthDistanceAndBearing[latb, longb, latc, longc]
dca = earthDistance[latc, longc, lata, longa]
println["dab: $dab"]
println["dbc: $dbc"]
println["dca: $dca"]
// Find bearings between seismic stations.
println["Bearing from a to b: " + (bab->DM)]
bac = earthBearing[lata, longa, latc, longc]
println["Bearing from a to c: " + (bac->DM)]
println["Bearing from b to c: " + (bbc->DM)]
// Times of arrival
ta1 = .1673 s
tb1 = .2890 s
tc1 = .5507 s
ta2 = .4240 s
tb2 = .6065 s
tc2 = .9990 s
// Velocities of waves
v1 = 6.0 km/s
v2 = 4.0 km/s
// This function was one I derived for finding the distance between
// the meteorite and the seismic station, given times of arrival for
// the P-waves and S-waves and their velocities.
d[t1, t2, v1=v1, v2=v2] := -((t1-t2) v1 v2)/(v1-v2)
// Find distances to seismic stations from the meteorite.
// I call the meteorite "point D"
dad = d[ta1, ta2]
dbd = d[tb1, tb2]
dcd = d[tc1, tc2]
println["dad: $dad"]
println["dbd: $dbd"]
println["dcd: $dcd"]
// Cosine rule for triangles... find the angle DAB between the meteor,
// point A and point B.
adab = arccos[(dad^2 + dab^2 - dbd^2)/(2 dad dab)]
println["Angle DAB: " + (adab -> DM)]
// You can either add or subtract this angle to the angle bab (the bearing
// from point a to point b.) I tried both ways, and only subtracting gave a
// self-consistent solution when all points were considered.
bad = bab - adab
println["Bearing from A to D: " + (bad -> DM)]
// Use calculation to find the resultant lat/long of the meteorite given
// the initial point (that of seismic station a), the distance, and the
// bearing of the meteorite.
[latd, longd] = resultantLatLong[lata, longa, dad, bad]
println["D is at lat: " + (latd->DM)]
println["D is at long: " + (longd->DM)]
// Reverse calculations. These distances and bearings should match the
// ones calculated above.
println["\nVerifying:"]
[dad, bad] = earthDistanceAndBearing[lata, longa, latd, longd]
println["Bearing from A to D: " + (bad -> DM)]
println["Distance from A to D: " + dad]
[dbd, bbd] = earthDistanceAndBearing[latb, longb, latd, longd]
println["Bearing from B to D: " + (bbd -> DM)]
println["Distance from B to D: " + dbd]
[dcd, bcd] = earthDistanceAndBearing[latc, longc, latd, longd]
println["Bearing from C to D: " + (bcd -> DM)]
println["Distance from C to D: " + dcd]
// Karen's solutions (assuming UTM zone 13, easting 519926, northing 4372660)
latk = DMS[39,30.199999] North
longk = DMS[104,46.093333] West
[diffDist, diffBearing] = earthDistanceAndBearing[latd, longd, latk, longk]
println["Difference in solutions is " + (diffDist->"ft")]
println["Bearing from Alan's solutions to Karen's: " + (diffBearing->DM)]
Download or view forensicSeismology.frink in plain text format
This is a program written in the programming language Frink.
For more information, view the Frink
Documentation or see More Sample Frink Programs.
Alan Eliasen was born 19966 days, 3 hours, 36 minutes ago.