SUBROUTINE GPT (DMJD,DLAT,DLON,DHGT,PRES,TEMP,UNDU) *+ * - - - - - - - - - * G P T * - - - - - - - - - * * This routine is part of the International Earth Rotation and * Reference Systems Service (IERS) Conventions software collection. * * This subroutine determines Global Pressure and Temperature (Boehm et al. 2007) * based on Spherical Harmonics up to degree and order 9. * * In general, Class 1, 2, and 3 models represent physical effects that * act on geodetic parameters while canonical models provide lower-level * representations or basic computations that are used by Class 1, 2, or * 3 models. * * Status: Class 1 model * * Class 1 models are those recommended to be used a priori in the * reduction of raw space geodetic data in order to determine * geodetic parameter estimates. * Class 2 models are those that eliminate an observational * singularity and are purely conventional in nature. * Class 3 models are those that are not required as either Class * 1 or 2. * Canonical models are accepted as is and cannot be classified as a * Class 1, 2, or 3 model. * * Given: * DMJD d Modified Julian Date * DLAT d Latitude given in radians (North Latitude) * DLON d Longitude given in radians (East Longitude) * DHGT d Ellipsoidal height in meters * * Returned: * PRES d Pressure given in hPa * TEMP d Temperature in degrees Celsius * UNDU d Geoid undulation in meters (Note 1) * * Notes: * * 1) This is from a 9x9 Earth Gravitational Model (EGM). * * Test case: * given input: DMJD = 55055D0 * DLAT = 0.6708665767D0 radians (NRAO, Green Bank, WV) * DLON = -1.393397187D0 radians * DHGT = 812.546 meters * expected output: PRES = 918.0710638757363995D0 hPa * TEMP = 19.31914181012882992D0 degrees Celsius * UNDU = -42.19185643717770517D0 meters * * References: * * Boehm, J., Heinkelmann, R. and Schuh, H., 2007, "Short Note: A * Global model of pressure and temperature for geodetic applications", * Journal of Geodesy, 81(10), pp. 679-683. * * Petit, G. and Luzum, B. (eds.), IERS Conventions (2010), * IERS Technical Note No. 36, BKG (2010) * * Revisions: * 2006 June 12 J. Boehm Original code * 2006 June 16 J. Boehm Accounted for geoid undulation * 2006 August 14 O. Montenbruck Recursions for Legendre polynomials * 2009 February 13 B.E. Stetzler Added header and copyright * 2009 March 30 B.E. Stetzler More modifications and defined twopi * 2009 March 31 B.E. Stetzler Provided test case * 2009 July 29 B.E. Stetzler Capitalized all variables for FORTRAN 77 * compatibility and corrected test case * latitude and longitude coordinates * 2011 October 18 B.E. Stetzler Corrected test case output mentioned * in the header *----------------------------------------------------------------------- IMPLICIT NONE DOUBLE PRECISION DMJD,DLAT,DLON,DHGT,PRES,TEMP,UNDU DOUBLE PRECISION V(10,10),W(10,10), . AP_MEAN(55),BP_MEAN(55),AP_AMP(55),BP_AMP(55), . AT_MEAN(55),BT_MEAN(55),AT_AMP(55),BT_AMP(55), . A_GEOID(55),B_GEOID(55) INTEGER*4 I,N,M,NMAX,MMAX DOUBLE PRECISION TWOPI,DOY,TEMP0,PRES0, . APM,APA,ATM,ATA,HORT,X,Y,Z PARAMETER (TWOPI = 6.283185307179586476925287D0) DATA (A_GEOID(I),I=1,55)/ .-5.6195D-001,-6.0794D-002,-2.0125D-001,-6.4180D-002,-3.6997D-002, .+1.0098D+001,+1.6436D+001,+1.4065D+001,+1.9881D+000,+6.4414D-001, .-4.7482D+000,-3.2290D+000,+5.0652D-001,+3.8279D-001,-2.6646D-002, .+1.7224D+000,-2.7970D-001,+6.8177D-001,-9.6658D-002,-1.5113D-002, .+2.9206D-003,-3.4621D+000,-3.8198D-001,+3.2306D-002,+6.9915D-003, .-2.3068D-003,-1.3548D-003,+4.7324D-006,+2.3527D+000,+1.2985D+000, .+2.1232D-001,+2.2571D-002,-3.7855D-003,+2.9449D-005,-1.6265D-004, .+1.1711D-007,+1.6732D+000,+1.9858D-001,+2.3975D-002,-9.0013D-004, .-2.2475D-003,-3.3095D-005,-1.2040D-005,+2.2010D-006,-1.0083D-006, .+8.6297D-001,+5.8231D-001,+2.0545D-002,-7.8110D-003,-1.4085D-004, .-8.8459D-006,+5.7256D-006,-1.5068D-006,+4.0095D-007,-2.4185D-008/ DATA (B_GEOID(I),I=1,55)/ .+0.0000D+000,+0.0000D+000,-6.5993D-002,+0.0000D+000,+6.5364D-002, .-5.8320D+000,+0.0000D+000,+1.6961D+000,-1.3557D+000,+1.2694D+000, .+0.0000D+000,-2.9310D+000,+9.4805D-001,-7.6243D-002,+4.1076D-002, .+0.0000D+000,-5.1808D-001,-3.4583D-001,-4.3632D-002,+2.2101D-003, .-1.0663D-002,+0.0000D+000,+1.0927D-001,-2.9463D-001,+1.4371D-003, .-1.1452D-002,-2.8156D-003,-3.5330D-004,+0.0000D+000,+4.4049D-001, .+5.5653D-002,-2.0396D-002,-1.7312D-003,+3.5805D-005,+7.2682D-005, .+2.2535D-006,+0.0000D+000,+1.9502D-002,+2.7919D-002,-8.1812D-003, .+4.4540D-004,+8.8663D-005,+5.5596D-005,+2.4826D-006,+1.0279D-006, .+0.0000D+000,+6.0529D-002,-3.5824D-002,-5.1367D-003,+3.0119D-005, .-2.9911D-005,+1.9844D-005,-1.2349D-006,-7.6756D-009,+5.0100D-008/ DATA (AP_MEAN(I),I=1,55)/ .+1.0108D+003,+8.4886D+000,+1.4799D+000,-1.3897D+001,+3.7516D-003, .-1.4936D-001,+1.2232D+001,-7.6615D-001,-6.7699D-002,+8.1002D-003, .-1.5874D+001,+3.6614D-001,-6.7807D-002,-3.6309D-003,+5.9966D-004, .+4.8163D+000,-3.7363D-001,-7.2071D-002,+1.9998D-003,-6.2385D-004, .-3.7916D-004,+4.7609D+000,-3.9534D-001,+8.6667D-003,+1.1569D-002, .+1.1441D-003,-1.4193D-004,-8.5723D-005,+6.5008D-001,-5.0889D-001, .-1.5754D-002,-2.8305D-003,+5.7458D-004,+3.2577D-005,-9.6052D-006, .-2.7974D-006,+1.3530D+000,-2.7271D-001,-3.0276D-004,+3.6286D-003, .-2.0398D-004,+1.5846D-005,-7.7787D-006,+1.1210D-006,+9.9020D-008, .+5.5046D-001,-2.7312D-001,+3.2532D-003,-2.4277D-003,+1.1596D-004, .+2.6421D-007,-1.3263D-006,+2.7322D-007,+1.4058D-007,+4.9414D-009/ DATA (BP_MEAN(I),I=1,55)/ .+0.0000D+000,+0.0000D+000,-1.2878D+000,+0.0000D+000,+7.0444D-001, .+3.3222D-001,+0.0000D+000,-2.9636D-001,+7.2248D-003,+7.9655D-003, .+0.0000D+000,+1.0854D+000,+1.1145D-002,-3.6513D-002,+3.1527D-003, .+0.0000D+000,-4.8434D-001,+5.2023D-002,-1.3091D-002,+1.8515D-003, .+1.5422D-004,+0.0000D+000,+6.8298D-001,+2.5261D-003,-9.9703D-004, .-1.0829D-003,+1.7688D-004,-3.1418D-005,+0.0000D+000,-3.7018D-001, .+4.3234D-002,+7.2559D-003,+3.1516D-004,+2.0024D-005,-8.0581D-006, .-2.3653D-006,+0.0000D+000,+1.0298D-001,-1.5086D-002,+5.6186D-003, .+3.2613D-005,+4.0567D-005,-1.3925D-006,-3.6219D-007,-2.0176D-008, .+0.0000D+000,-1.8364D-001,+1.8508D-002,+7.5016D-004,-9.6139D-005, .-3.1995D-006,+1.3868D-007,-1.9486D-007,+3.0165D-010,-6.4376D-010/ DATA (AP_AMP(I),I=1,55)/ .-1.0444D-001,+1.6618D-001,-6.3974D-002,+1.0922D+000,+5.7472D-001, .-3.0277D-001,-3.5087D+000,+7.1264D-003,-1.4030D-001,+3.7050D-002, .+4.0208D-001,-3.0431D-001,-1.3292D-001,+4.6746D-003,-1.5902D-004, .+2.8624D+000,-3.9315D-001,-6.4371D-002,+1.6444D-002,-2.3403D-003, .+4.2127D-005,+1.9945D+000,-6.0907D-001,-3.5386D-002,-1.0910D-003, .-1.2799D-004,+4.0970D-005,+2.2131D-005,-5.3292D-001,-2.9765D-001, .-3.2877D-002,+1.7691D-003,+5.9692D-005,+3.1725D-005,+2.0741D-005, .-3.7622D-007,+2.6372D+000,-3.1165D-001,+1.6439D-002,+2.1633D-004, .+1.7485D-004,+2.1587D-005,+6.1064D-006,-1.3755D-008,-7.8748D-008, .-5.9152D-001,-1.7676D-001,+8.1807D-003,+1.0445D-003,+2.3432D-004, .+9.3421D-006,+2.8104D-006,-1.5788D-007,-3.0648D-008,+2.6421D-010/ DATA (BP_AMP(I),I=1,55)/ .+0.0000D+000,+0.0000D+000,+9.3340D-001,+0.0000D+000,+8.2346D-001, .+2.2082D-001,+0.0000D+000,+9.6177D-001,-1.5650D-002,+1.2708D-003, .+0.0000D+000,-3.9913D-001,+2.8020D-002,+2.8334D-002,+8.5980D-004, .+0.0000D+000,+3.0545D-001,-2.1691D-002,+6.4067D-004,-3.6528D-005, .-1.1166D-004,+0.0000D+000,-7.6974D-002,-1.8986D-002,+5.6896D-003, .-2.4159D-004,-2.3033D-004,-9.6783D-006,+0.0000D+000,-1.0218D-001, .-1.3916D-002,-4.1025D-003,-5.1340D-005,-7.0114D-005,-3.3152D-007, .+1.6901D-006,+0.0000D+000,-1.2422D-002,+2.5072D-003,+1.1205D-003, .-1.3034D-004,-2.3971D-005,-2.6622D-006,+5.7852D-007,+4.5847D-008, .+0.0000D+000,+4.4777D-002,-3.0421D-003,+2.6062D-005,-7.2421D-005, .+1.9119D-006,+3.9236D-007,+2.2390D-007,+2.9765D-009,-4.6452D-009/ DATA (AT_MEAN(I),I=1,55)/ .+1.6257D+001,+2.1224D+000,+9.2569D-001,-2.5974D+001,+1.4510D+000, .+9.2468D-002,-5.3192D-001,+2.1094D-001,-6.9210D-002,-3.4060D-002, .-4.6569D+000,+2.6385D-001,-3.6093D-002,+1.0198D-002,-1.8783D-003, .+7.4983D-001,+1.1741D-001,+3.9940D-002,+5.1348D-003,+5.9111D-003, .+8.6133D-006,+6.3057D-001,+1.5203D-001,+3.9702D-002,+4.6334D-003, .+2.4406D-004,+1.5189D-004,+1.9581D-007,+5.4414D-001,+3.5722D-001, .+5.2763D-002,+4.1147D-003,-2.7239D-004,-5.9957D-005,+1.6394D-006, .-7.3045D-007,-2.9394D+000,+5.5579D-002,+1.8852D-002,+3.4272D-003, .-2.3193D-005,-2.9349D-005,+3.6397D-007,+2.0490D-006,-6.4719D-008, .-5.2225D-001,+2.0799D-001,+1.3477D-003,+3.1613D-004,-2.2285D-004, .-1.8137D-005,-1.5177D-007,+6.1343D-007,+7.8566D-008,+1.0749D-009/ DATA (BT_MEAN(I),I=1,55)/ .+0.0000D+000,+0.0000D+000,+1.0210D+000,+0.0000D+000,+6.0194D-001, .+1.2292D-001,+0.0000D+000,-4.2184D-001,+1.8230D-001,+4.2329D-002, .+0.0000D+000,+9.3312D-002,+9.5346D-002,-1.9724D-003,+5.8776D-003, .+0.0000D+000,-2.0940D-001,+3.4199D-002,-5.7672D-003,-2.1590D-003, .+5.6815D-004,+0.0000D+000,+2.2858D-001,+1.2283D-002,-9.3679D-003, .-1.4233D-003,-1.5962D-004,+4.0160D-005,+0.0000D+000,+3.6353D-002, .-9.4263D-004,-3.6762D-003,+5.8608D-005,-2.6391D-005,+3.2095D-006, .-1.1605D-006,+0.0000D+000,+1.6306D-001,+1.3293D-002,-1.1395D-003, .+5.1097D-005,+3.3977D-005,+7.6449D-006,-1.7602D-007,-7.6558D-008, .+0.0000D+000,-4.5415D-002,-1.8027D-002,+3.6561D-004,-1.1274D-004, .+1.3047D-005,+2.0001D-006,-1.5152D-007,-2.7807D-008,+7.7491D-009/ DATA (AT_AMP(I),I=1,55)/ .-1.8654D+000,-9.0041D+000,-1.2974D-001,-3.6053D+000,+2.0284D-002, .+2.1872D-001,-1.3015D+000,+4.0355D-001,+2.2216D-001,-4.0605D-003, .+1.9623D+000,+4.2887D-001,+2.1437D-001,-1.0061D-002,-1.1368D-003, .-6.9235D-002,+5.6758D-001,+1.1917D-001,-7.0765D-003,+3.0017D-004, .+3.0601D-004,+1.6559D+000,+2.0722D-001,+6.0013D-002,+1.7023D-004, .-9.2424D-004,+1.1269D-005,-6.9911D-006,-2.0886D+000,-6.7879D-002, .-8.5922D-004,-1.6087D-003,-4.5549D-005,+3.3178D-005,-6.1715D-006, .-1.4446D-006,-3.7210D-001,+1.5775D-001,-1.7827D-003,-4.4396D-004, .+2.2844D-004,-1.1215D-005,-2.1120D-006,-9.6421D-007,-1.4170D-008, .+7.8720D-001,-4.4238D-002,-1.5120D-003,-9.4119D-004,+4.0645D-006, .-4.9253D-006,-1.8656D-006,-4.0736D-007,-4.9594D-008,+1.6134D-009/ DATA (BT_AMP(I),I=1,55)/ .+0.0000D+000,+0.0000D+000,-8.9895D-001,+0.0000D+000,-1.0790D+000, .-1.2699D-001,+0.0000D+000,-5.9033D-001,+3.4865D-002,-3.2614D-002, .+0.0000D+000,-2.4310D-002,+1.5607D-002,-2.9833D-002,-5.9048D-003, .+0.0000D+000,+2.8383D-001,+4.0509D-002,-1.8834D-002,-1.2654D-003, .-1.3794D-004,+0.0000D+000,+1.3306D-001,+3.4960D-002,-3.6799D-003, .-3.5626D-004,+1.4814D-004,+3.7932D-006,+0.0000D+000,+2.0801D-001, .+6.5640D-003,-3.4893D-003,-2.7395D-004,+7.4296D-005,-7.9927D-006, .-1.0277D-006,+0.0000D+000,+3.6515D-002,-7.4319D-003,-6.2873D-004, .-8.2461D-005,+3.1095D-005,-5.3860D-007,-1.2055D-007,-1.1517D-007, .+0.0000D+000,+3.1404D-002,+1.5580D-002,-1.1428D-003,+3.3529D-005, .+1.0387D-005,-1.9378D-006,-2.7327D-007,+7.5833D-009,-9.2323D-009/ * Reference day is 28 January 1980 * This is taken from Niell (1996) to be consistent (See References) * For constant values use: doy = 91.3125 DOY = DMJD - 44239D0 + 1 - 28 * Define degree n and order m EGM NMAX = 9 MMAX = 9 * Define unit vector X = DCOS(DLAT)*DCOS(DLON) Y = DCOS(DLAT)*DSIN(DLON) Z = DSIN(DLAT) * Legendre polynomials V(1,1) = 1.0D0 W(1,1) = 0.0D0; V(2,1) = Z * V(1,1); W(2,1) = 0.0; DO N=2,NMAX V(N+1,1) = ((2*N-1) * Z * V(N,1) - (N-1) * V(N-1,1)) / N W(N+1,1) = 0.0D0 ENDDO DO M=1,NMAX V(M+1,M+1) = (2*M-1) * (X*V(M,M) - Y*W(M,M)) W(M+1,M+1) = (2*M-1) * (X*W(M,M) + Y*V(M,M)) IF (M < NMAX) THEN V(M+2,M+1) = (2*M+1) * Z * V(M+1,M+1) W(M+2,M+1) = (2*M+1) * Z * W(M+1,M+1) ENDIF DO N=M+2,NMAX V(N+1,M+1) = ((2*N-1)*Z*V(N,M+1) - (N+M-1)*V(N-1,M+1)) / (N-M) W(N+1,M+1) = ((2*N-1)*Z*W(N,M+1) - (N+M-1)*W(N-1,M+1)) / (N-M) ENDDO ENDDO * Geoidal height UNDU = 0.D0 I = 0 DO N=0,NMAX DO M=0,N I = I+1 UNDU = UNDU + (A_GEOID(I)*V(N+1,M+1) + B_GEOID(I)*W(N+1,M+1)) ENDDO ENDDO * orthometric height HORT = DHGT - UNDU * Surface pressure on the geoid APM = 0.D0 APA = 0.D0 I = 0 DO N=0,NMAX DO M=0,N I = I+1 APM = APM + (AP_MEAN(I)*V(N+1,M+1) + BP_MEAN(I)*W(N+1,M+1)) APA = APA + (AP_AMP(I) *V(N+1,M+1) + BP_AMP(I) *W(N+1,M+1)) ENDDO ENDDO PRES0 = APM + APA*DCOS(DOY/365.25D0*TWOPI) * height correction for pressure PRES = PRES0*(1D0-0.0000226D0*HORT)**5.225D0 * Surface temperature on the geoid ATM = 0D0 ATA = 0D0 I = 0 DO N=0,NMAX DO M=0,N I = I+1 ATM = ATM + (AT_MEAN(I)*V(N+1,M+1) + BT_MEAN(I)*W(N+1,M+1)) ATA = ATA + (AT_AMP(I) *V(N+1,M+1) + BT_AMP(I) *W(N+1,M+1)) ENDDO ENDDO TEMP0 = ATM + ATA*DCOS(DOY/365.25D0*TWOPI) * height correction for temperature TEMP = TEMP0 - 0.0065D0*HORT * Finished. *+---------------------------------------------------------------------- * * Copyright (C) 2008 * IERS Conventions Center * * ================================== * IERS Conventions Software License * ================================== * * NOTICE TO USER: * * BY USING THIS SOFTWARE YOU ACCEPT THE FOLLOWING TERMS AND CONDITIONS * WHICH APPLY TO ITS USE. * * 1. The Software is provided by the IERS Conventions Center ("the * Center"). * * 2. Permission is granted to anyone to use the Software for any * purpose, including commercial applications, free of charge, * subject to the conditions and restrictions listed below. * * 3. You (the user) may adapt the Software and its algorithms for your * own purposes and you may distribute the resulting "derived work" * to others, provided that the derived work complies with the * following requirements: * * a) Your work shall be clearly identified so that it cannot be * mistaken for IERS Conventions software and that it has been * neither distributed by nor endorsed by the Center. * * b) Your work (including source code) must contain descriptions of * how the derived work is based upon and/or differs from the * original Software. * * c) The name(s) of all modified routine(s) that you distribute * shall be changed. * * d) The origin of the IERS Conventions components of your derived * work must not be misrepresented; you must not claim that you * wrote the original Software. * * e) The source code must be included for all routine(s) that you * distribute. This notice must be reproduced intact in any * source distribution. * * 4. In any published work produced by the user and which includes * results achieved by using the Software, you shall acknowledge * that the Software was used in obtaining those results. * * 5. The Software is provided to the user "as is" and the Center makes * no warranty as to its use or performance. The Center does not * and cannot warrant the performance or results which the user may * obtain by using the Software. The Center makes no warranties, * express or implied, as to non-infringement of third party rights, * merchantability, or fitness for any particular purpose. In no * event will the Center be liable to the user for any consequential, * incidental, or special damages, including any lost profits or lost * savings, even if a Center representative has been advised of such * damages, or for any claim by any third party. * * Correspondence concerning IERS Conventions software should be * addressed as follows: * * Gerard Petit * Internet email: gpetit[at]bipm.org * Postal address: IERS Conventions Center * Time, frequency and gravimetry section, BIPM * Pavillon de Breteuil * 92312 Sevres FRANCE * * or * * Brian Luzum * Internet email: brian.luzum[at]usno.navy.mil * Postal address: IERS Conventions Center * Earth Orientation Department * 3450 Massachusetts Ave, NW * Washington, DC 20392 * * *----------------------------------------------------------------------- END