|
SURFEX v7.3
General documentation of Surfex
|
|
SURFEX v7.3
General documentation of Surfex
|
00001 ! ################################################### 00002 SUBROUTINE OI_LATLON_CONF_PROJ(NDIM,PLAT0,PLON0,PRPK,PBETA,PLATOR,PLONOR, & 00003 PX,PY,PLAT,PLON ) 00004 ! ################################################### 00005 ! 00006 !!**** *LATLON_CONF_PROJ * - Routine to compute geographical coordinates 00007 !! 00008 !! PURPOSE 00009 !! ------- 00010 ! This routine computes the latitude and longitude of 00011 ! an array given in cartesian conformal coordinates 00012 ! Five map projections are available: 00013 ! - polar-stereographic from south pole (PRPK=1), 00014 ! - lambert conformal from south pole (0<PRPK<1), 00015 ! - mercator (PRPK=0), 00016 ! - lambert conformal from north pole (-1<PRPK<0), 00017 ! - polar-stereographic from north pole (PRPK=-1). 00018 ! 00019 ! 00020 !!** METHOD 00021 !! ------ 00022 !! Spherical earth approximation is used. Longitude origin is 00023 !! set in Greenwich, and is positive eastwards. An anticlockwise 00024 !! rotation of PBETA degrees is applied to the conformal frame 00025 !! with respect to the geographical directions. 00026 !! 00027 !! WARNING: ALL INPUT AND OUTPUT ANGLES ARE IN DEGREES... 00028 !! 00029 !! EXTERNAL 00030 !! -------- 00031 !! None 00032 !! 00033 !! REFERENCE 00034 !! --------- 00035 !! Asencio N. et al., 1994, "Le projet de modele non-hydrostatique 00036 !! commun CNRM-LA, specifications techniques", 00037 !! Note CNRM/GMME, 26, 139p, (Chapter 2). 00038 !! Ducrocq V., 1994, "Generation de la grille dans le modele", 00039 !! Note interne MNH, 5 mai, 3p. 00040 !! Joly A., 1992, "Geographic parameters for ARPEGE/ALADIN", 00041 !! Internal note ARPEGE/ALADIN, february 27,28p. 00042 !! Levallois J., 1970, "Geodesie generale", Tome 2, Collection 00043 !! de l'IGN, Eyrolles, Paris, 408p. 00044 !! 00045 !! AUTHOR 00046 !! ------ 00047 !! P.M. *LA* 00048 !! 00049 !! MODIFICATION 00050 !! ------------ 00051 !! Original PM 24/05/94 00052 !! Updated PM 27/07/94 00053 !! Updated VD 23/08/94 00054 !! Updated VM 24/10/95 projection from north pole (PRPK<0) and 00055 !! longitudes set between PLON0-180. and PLON0+180. 00056 !! Update VM 11/2004 externalized version 00057 !------------------------------------------------------------------------------- 00058 ! 00059 !* 0. DECLARATIONS 00060 ! ------------ 00061 ! 00062 USE MODD_CSTS,ONLY : XPI, XRADIUS 00063 ! 00064 ! 00065 USE YOMHOOK ,ONLY : LHOOK, DR_HOOK 00066 USE PARKIND1 ,ONLY : JPRB 00067 ! 00068 IMPLICIT NONE 00069 ! 00070 !* 0.1 Declarations of arguments and results 00071 ! 00072 REAL, INTENT(IN) :: PLAT0 ! Reference latitude 00073 REAL, INTENT(IN) :: PLON0 ! Reference longitude 00074 REAL, INTENT(IN) :: PRPK ! projection parameter 00075 ! ! K=1 : stereographic north pole 00076 ! ! 0<K<1 : Lambert, north hemisphere 00077 ! ! K=0 : Mercator 00078 ! !-1<K<0 : Lambert, south hemisphere 00079 ! ! K=-1: stereographic south pole 00080 REAL, INTENT(IN) :: PBETA ! angle between grid and reference longitude 00081 REAL, INTENT(IN) :: PLATOR ! Latitude of the origine point 00082 REAL, INTENT(IN) :: PLONOR ! Longitude of the origine point 00083 REAL, DIMENSION(NDIM), INTENT(IN) :: PX,PY 00084 ! given conformal coordinates of the 00085 ! processed points (meters); 00086 REAL, DIMENSION(NDIM), INTENT(OUT):: PLAT,PLON 00087 ! returned geographic latitudes and 00088 ! longitudes of the processed points 00089 ! (degrees). 00090 INTEGER, INTENT(IN) :: NDIM 00091 ! 00092 !* 0.2 Declarations of local variables 00093 ! 00094 REAL, DIMENSION(NDIM) :: ZY 00095 REAL :: ZRPK,ZBETA,ZLAT0,ZLON0,ZLATOR,ZLONOR 00096 REAL :: ZRDSDG,ZCLAT0,ZSLAT0,ZCLATOR,ZSLATOR 00097 REAL :: ZXBM0,ZYBM0,ZRO0,ZGA0 00098 REAL :: ZXP,ZYP,ZEPSI,ZT1,ZCGAM,ZSGAM,ZRACLAT0 00099 ! 00100 REAL, DIMENSION(NDIM) :: ZATA,ZRO2,ZT2,ZXMI0,ZYMI0 00101 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00102 ! 00103 !-------------------------------------------------------------------------------- 00104 ! 00105 !* 1. Preliminary calculations for all projections 00106 ! -------------------------------------------- 00107 ! 00108 IF (LHOOK) CALL DR_HOOK('OI_LATLON_CONF_PROJ',0,ZHOOK_HANDLE) 00109 ZRDSDG = XPI/180. ! Degree to radian conversion factor 00110 ZEPSI = 10.*EPSILON(1.) ! A small number 00111 ! 00112 ! By definition, (PLONOR,PLATOR) are the geographical 00113 ! coordinates, and (ZXBM0,ZYBM0) the conformal cartesian 00114 ! point of coordinates (x=0,y=0). 00115 ! 00116 ZXBM0 = 0. 00117 ZYBM0 = 0. 00118 ! 00119 !------------------------------------------------------------------------------- 00120 ! 00121 !* 2. POLAR STEREOGRAPHIC AND LAMBERT CONFORMAL CASES 00122 ! ----------------------------------------------- 00123 ! (PRPK=1 P-stereo, 0<PRPK<1 Lambert) 00124 ! 00125 IF(PRPK /= 0.) THEN 00126 ! 00127 IF (PRPK<0.) THEN ! projection from north pole 00128 ZRPK=-PRPK 00129 ZBETA=-PBETA 00130 ZLAT0=-PLAT0 00131 ZLON0=PLON0+180. 00132 ZLATOR=-PLATOR 00133 ZLONOR=PLONOR+180. 00134 ZY(:)=-PY(:) 00135 ZYBM0=-ZYBM0 00136 ELSE ! projection from south pole 00137 ZRPK=PRPK 00138 ZBETA=PBETA 00139 ZLAT0=PLAT0 00140 ZLON0=PLON0 00141 ZLATOR=PLATOR 00142 ZLONOR=PLONOR 00143 ZY(:)=PY(:) 00144 ENDIF 00145 ! 00146 !* 2.1 Preliminary calculations 00147 ! 00148 ZCLAT0 = COS(ZRDSDG*ZLAT0) 00149 ZSLAT0 = SIN(ZRDSDG*ZLAT0) 00150 ZCLATOR = COS(ZRDSDG*ZLATOR) 00151 ZSLATOR = SIN(ZRDSDG*ZLATOR) 00152 ZRO0 = (XRADIUS/ZRPK)*(ABS(ZCLAT0))**(1.-ZRPK) & 00153 * ((1.+ZSLAT0)*ABS(ZCLATOR)/(1.+ZSLATOR))**ZRPK 00154 ZGA0 = (ZRPK*(ZLONOR-ZLON0)-ZBETA)*ZRDSDG 00155 ZXP = ZXBM0-ZRO0*SIN(ZGA0) 00156 ZYP = ZYBM0+ZRO0*COS(ZGA0) 00157 ! 00158 !* 2.2 Longitude 00159 ! 00160 WHERE (ABS(ZY(:)-ZYP) < ZEPSI & 00161 .AND.ABS(PX(:)-ZXP) < ZEPSI) 00162 ZATA(:) = 0. 00163 ELSEWHERE 00164 ZATA(:) = ATAN2(-(ZXP-PX(:)),(ZYP-ZY(:)))/ZRDSDG 00165 END WHERE 00166 ! 00167 PLON(:) = (ZBETA+ZATA(:))/ZRPK+ZLON0 00168 ! 00169 !* 2.3 Latitude 00170 ! 00171 ZRO2(:) = (PX(:)-ZXP)**2+(ZY(:)-ZYP)**2 00172 ZT1 = (XRADIUS*(ABS(ZCLAT0))**(1.-ZRPK))**(2./ZRPK) & 00173 * (1+ZSLAT0)**2 00174 ZT2(:) = (ZRPK**2*ZRO2(:))**(1./ZRPK) 00175 ! 00176 PLAT(:) = (XPI/2.-ACOS((ZT1-ZT2(:))/(ZT1+ZT2(:))))/ZRDSDG 00177 ! 00178 IF (PRPK<0.) THEN ! projection from north pole 00179 PLAT(:)=-PLAT(:) 00180 PLON(:)=PLON(:)+180. 00181 ENDIF 00182 ! 00183 !------------------------------------------------------------------------------- 00184 ! 00185 !* 3. MERCATOR PROJECTION WITH ROTATION 00186 ! --------------------------------- 00187 ! (PRPK=0) 00188 ! 00189 ELSE 00190 ! 00191 !* 3.1 Preliminary calculations 00192 ! 00193 ZCGAM = COS(-ZRDSDG*PBETA) 00194 ZSGAM = SIN(-ZRDSDG*PBETA) 00195 ZRACLAT0 = XRADIUS*COS(ZRDSDG*PLAT0) 00196 ! 00197 !* 3.2 Longitude 00198 ! 00199 ZXMI0(:) = PX(:)-ZXBM0 00200 ZYMI0(:) = PY(:)-ZYBM0 00201 ! 00202 PLON(:) = (ZXMI0(:)*ZCGAM+ZYMI0(:)*ZSGAM) & 00203 / (ZRACLAT0*ZRDSDG)+PLONOR 00204 ! 00205 !* 3.3 Latitude 00206 ! 00207 ZT1 = ALOG(TAN(XPI/4.+PLATOR*ZRDSDG/2.)) 00208 ZT2(:) = (-ZXMI0(:)*ZSGAM+ZYMI0(:)*ZCGAM)/ZRACLAT0 00209 ! 00210 PLAT(:) = (-XPI/2.+2.*ATAN(EXP(ZT1+ZT2(:))))/ZRDSDG 00211 ! 00212 !--------------------------------------------------------------------------------- 00213 ! 00214 !* 4. EXIT 00215 ! ---- 00216 ! 00217 END IF 00218 PLON(:)=PLON(:)+NINT((PLON0-PLON(:))/360.)*360. 00219 IF (LHOOK) CALL DR_HOOK('OI_LATLON_CONF_PROJ',1,ZHOOK_HANDLE) 00220 !--------------------------------------------------------------------------------- 00221 END SUBROUTINE OI_LATLON_CONF_PROJ
1.8.0
