W3FB10TheW3FB10routine computes the bearing and great circle distance. from point (1) to point (2) assuming a spherical earth. The North and South poles are special cases. If latitude of point (1) is within 1e-10 degrees of the North pole, bearing is the negative longitude of point (2) by convention. If latitude of point (1) is within 1e-10 degrees of the south pole, bearing is the longitude of point (2) by convention. If point (2) is within 1e-6 radians of the antipode of point (1), the bearing will be set to zero. If point (1) and point (2) are within 1e-10 radians of each other, both bearing and distance will be set to zero. USAGE: CALL W3FB10(DLAT1, DLON1, DLAT2, DLON2, BEARD, GCDKM) Input argument list: DLAT1 - REAL LATITUDE OF POINT (1) IN DEGREES NORTH. DLON1 - REAL LONGITUDE OF POINT (1) IN DEGREES EAST. DLAT2 - REAL LATITUDE OF POINT (2) IN DEGREES NORTH. DLON2 - REAL LONGITUDE OF POINT (2) IN DEGREES EAST. Output argument list: BEARD - REAL BEARING OF POINT (2) FROM POINT (1) IN COMPASS DEGREES WITH NORTH = 0.0, VALUES FROM -180.0 TO +180.0 DEGREES. GCDKM - REAL GREAT CIRCLE DISTANCE FROM POINT (1) TO POINT (2) IN KILOMETERS. REMARKS: According to the NMC Handbook, the earth's radius is 6371.2 kilometers. This is what we use, even though the value recommended by the Smithsonian Meteorological Handbook is 6371.221 KM. (I wouldn't want you to think that i didn't know what the correct value was.) METHOD: The poles are special cases, and handled separately. Otherwise, from spherical trigonometry, the law of cosines is used to calculate the third side of the spherical triangle having sides from the pole to points (1) and (2) (the colatitudes). Then the law of sines is used to calculate the angle at point (1). A test is applied to see whether the arcsine result may be be used as such, giving an acute angle as the bearing, or whether The arcsine result should be subtracted from pi, giving an obtuse angle as the bearing. This test is derived by constructing a right spherical triangle using the pole, point (2), and the meridian through point(1). The latitude of the right-angled vertex then provides a test--if latitude (1) is greater than this latitude, the bearing angle must be obtuse, otherwise acute. If the two points are within 1e-6 radians of each other a flat earth is assumed, and the four-quadrant arctangent function is used to find the bearing. The y-displacement is the difference in latitude and the x-displacement is the difference in longitude times cosine latitude, both in radians. distance is then the diagonal. Fundamental trigonometric identities are used freely, such as that cos(x) = sin(pi/2 - x), etc. See almost any mathematical Handbook, such as the C.R.C. standard math tables under 'relations in any spherical triangle', or the National Bureau of Standards 'Handbook of mathematical functions' under section 4.3.149, formulas for solution of spherical triangles. Double precision is used internally because of the wide range of geographic values that may be used.

W3lib.tar Library contains Fortran 90 decoder/encoder routines for GRIB edition 1. (Fortran90)

Date posted: 2/22/2007