|
SURFEX v7.3
General documentation of Surfex
|
|
SURFEX v7.3
General documentation of Surfex
|
00001 ! ######### 00002 SUBROUTINE HOR_INTERPOL_GAUSS(KLUOUT,PFIELDIN,PFIELDOUT) 00003 ! ################################################################################# 00004 ! 00005 !!**** *HOR_INTERPOL_GAUSS* - Interpolation from a gaussian grid 00006 !! 00007 !! PURPOSE 00008 !! ------- 00009 ! 00010 !!** METHOD 00011 !! ------ 00012 !! 00013 !! REFERENCE 00014 !! --------- 00015 !! 00016 !! 00017 !! AUTHOR 00018 !! ------ 00019 !! V. Masson 00020 !! 00021 !! MODIFICATIONS 00022 !! ------------- 00023 !! Original 01/2004 00024 !!------------------------------------------------------------------ 00025 ! 00026 ! 00027 ! 00028 USE MODD_PREP, ONLY : XLAT_OUT, XLON_OUT, LINTERP 00029 USE MODD_GRID_GAUSS, ONLY : XILA1, XILO1, XILA2, XILO2, NINLA, NINLO, NILEN, LROTPOLE, & 00030 XLAP, XLOP, XCOEF 00031 USE MODD_GRID_GRIB, ONLY : NNI 00032 USE MODD_SURF_PAR, ONLY : XUNDEF 00033 ! 00034 USE MODI_HORIBL_SURF 00035 ! 00036 ! 00037 USE YOMHOOK ,ONLY : LHOOK, DR_HOOK 00038 USE PARKIND1 ,ONLY : JPRB 00039 ! 00040 IMPLICIT NONE 00041 ! 00042 !* 0.1 declarations of arguments 00043 ! 00044 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00045 REAL, DIMENSION(:,:), INTENT(IN) :: PFIELDIN ! field to interpolate horizontally 00046 REAL, DIMENSION(:,:), INTENT(OUT) :: PFIELDOUT ! interpolated field 00047 ! 00048 !* 0.2 declarations of local variables 00049 ! 00050 REAL, DIMENSION(:), ALLOCATABLE :: ZLAT ! latitudes 00051 REAL, DIMENSION(:), ALLOCATABLE :: ZLON ! longitudes 00052 INTEGER, DIMENSION(:), ALLOCATABLE :: IMASKIN ! input mask 00053 INTEGER, DIMENSION(:), ALLOCATABLE :: IMASKOUT ! output mask 00054 INTEGER :: INO ! output number of points 00055 INTEGER :: JL ! loop counter 00056 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00057 ! 00058 !------------------------------------------------------------------------------------- 00059 ! 00060 !* 1. Allocations 00061 ! 00062 IF (LHOOK) CALL DR_HOOK('HOR_INTERPOL_GAUSS',0,ZHOOK_HANDLE) 00063 INO = SIZE(XLAT_OUT) 00064 ! 00065 ALLOCATE(IMASKIN (NNI)) 00066 ! 00067 ALLOCATE(ZLAT (INO)) 00068 ALLOCATE(ZLON (INO)) 00069 ALLOCATE(IMASKOUT(INO)) 00070 IMASKOUT = 1 00071 ! 00072 !* 2. Is pole rotated? 00073 ! 00074 IF (LROTPOLE) THEN 00075 !* transformation of output latitudes, longitudes into rotated coordinates 00076 CALL ARPEGE_STRETCH_A(INO,XLAP,XLOP,XCOEF,XLAT_OUT,XLON_OUT,ZLAT,ZLON) 00077 ELSE 00078 ZLAT = XLAT_OUT 00079 ZLON = XLON_OUT 00080 END IF 00081 ! 00082 !* 3. Input mask 00083 ! 00084 DO JL=1,SIZE(PFIELDIN,2) 00085 ! 00086 IMASKIN(:) = 1. 00087 WHERE(PFIELDIN(:,JL)==XUNDEF) IMASKIN = 0. 00088 ! 00089 ! 00090 !* 4. Interpolation with horibl 00091 ! 00092 CALL HORIBL_SURF(XILA1,XILO1,XILA2,XILO2,NINLA,NINLO,NNI,PFIELDIN(:,JL),INO,ZLON,ZLAT,PFIELDOUT(:,JL), & 00093 .FALSE.,KLUOUT,LINTERP,IMASKIN,IMASKOUT) 00094 END DO 00095 ! 00096 ! 00097 !* 5. Deallocations 00098 ! 00099 DEALLOCATE(ZLAT) 00100 DEALLOCATE(ZLON) 00101 DEALLOCATE(IMASKIN ) 00102 DEALLOCATE(IMASKOUT) 00103 ! 00104 !------------------------------------------------------------------------------- 00105 !------------------------------------------------------------------------------- 00106 IF (LHOOK) CALL DR_HOOK('HOR_INTERPOL_GAUSS',1,ZHOOK_HANDLE) 00107 CONTAINS 00108 ! 00109 ! ########################################################################## 00110 SUBROUTINE ARPEGE_STRETCH_A(KN,PLAP,PLOP,PCOEF,PLAR,PLOR,PLAC,PLOC) 00111 ! ########################################################################## 00112 !!**** *ARPEGE_STRETCH_A* - Projection to Arpege stretched grid 00113 !! 00114 !! PURPOSE 00115 !! ------- 00116 !! 00117 !! Projection from standard Lat,Lon grid to Arpege stretched grid 00118 !! 00119 !! METHOD 00120 !! ------ 00121 !! 00122 !! The projection is defined in two steps : 00123 !! 1. A rotation to place the stretching pole at the north pole 00124 !! 2. The stretching 00125 !! This routine is a basic implementation of the informations founded in 00126 !! 'Note de travail Arpege n#3' 00127 !! 'Transformation de coordonnees' 00128 !! J.F.Geleyn 1988 00129 !! This document describes a slightly different transformation in 3 steps. Only the 00130 !! two first steps are to be taken in account (at the time of writing this paper has 00131 !! not been updated). 00132 !! 00133 !! EXTERNAL 00134 !! -------- 00135 !! 00136 !! 00137 !! IMPLICIT ARGUMENTS 00138 !! ------------------ 00139 !! 00140 !! REFERENCE 00141 !! --------- 00142 !! 00143 !! This routine is based on : 00144 !! 'Note de travail ARPEGE' number 3 00145 !! by J.F. GELEYN (may 1988) 00146 !! 00147 !! AUTHOR 00148 !! ------ 00149 !! 00150 !! V.Bousquet 00151 !! 00152 !! MODIFICATIONS 00153 !! ------------- 00154 !! 00155 !! Original 07/01/1999 00156 !! V. Masson 01/2004 Externalization of surface 00157 !! 00158 ! 00159 ! 0. DECLARATIONS 00160 ! --------------- 00161 ! 00162 USE MODD_CSTS,ONLY : XPI 00163 ! 00164 IMPLICIT NONE 00165 ! 00166 ! 0.1. Declaration of arguments 00167 ! ----------------------------- 00168 00169 INTEGER, INTENT(IN) :: KN ! Number of points to convert 00170 REAL, INTENT(IN) :: PLAP ! Latitude of stretching pole 00171 REAL, INTENT(IN) :: PLOP ! Longitude of stretching pole 00172 REAL, INTENT(IN) :: PCOEF ! Stretching coefficient 00173 REAL, DIMENSION(KN), INTENT(IN) :: PLAR ! Lat. of points 00174 REAL, DIMENSION(KN), INTENT(IN) :: PLOR ! Lon. of points 00175 REAL, DIMENSION(KN), INTENT(OUT) :: PLAC ! Computed pseudo-lat. of points 00176 REAL, DIMENSION(KN), INTENT(OUT) :: PLOC ! Computed pseudo-lon. of points 00177 ! 00178 REAL :: ZSINSTRETCHLA ! Sine of stretching point lat. 00179 REAL :: ZSINSTRETCHLO ! Sine of stretching point lon. 00180 REAL :: ZCOSSTRETCHLA ! Cosine of stretching point lat. 00181 REAL :: ZCOSSTRETCHLO ! Cosine of stretching point lon. 00182 REAL :: ZSINLA ! Sine of computed point latitude 00183 REAL :: ZSINLO ! Sine of computed point longitude 00184 REAL :: ZCOSLA ! Cosine of computed point latitude 00185 REAL :: ZCOSLO ! Cosine of computed point longitude 00186 REAL :: ZSINLAS ! Sine of point's pseudo-latitude 00187 REAL :: ZSINLOS ! Sine of point's pseudo-longitude 00188 REAL :: ZCOSLOS ! Cosine of point's pseudo-lon. 00189 REAL :: ZA,ZB,ZD ! Dummy variables used for 00190 REAL :: ZX,ZY ! computations 00191 ! 00192 INTEGER :: JLOOP1 ! Dummy loop counter 00193 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00194 ! 00195 IF (LHOOK) CALL DR_HOOK('ARPEGE_STRETCH_A',0,ZHOOK_HANDLE) 00196 ZSINSTRETCHLA = SIN(PLAP*XPI/180.) 00197 ZCOSSTRETCHLA = COS(PLAP*XPI/180.) 00198 ZSINSTRETCHLO = SIN(PLOP*XPI/180.) 00199 ZCOSSTRETCHLO = COS(PLOP*XPI/180.) 00200 ! L = longitude (0 = Greenwich, + toward east) 00201 ! l = latitude (90 = N.P., -90 = S.P.) 00202 ! p stands for stretching pole 00203 PLAC(:) = PLAR(:) * XPI / 180. 00204 PLOC(:) = PLOR(:) * XPI / 180. 00205 ! A = 1 + c.c 00206 ZA = 1. + PCOEF*PCOEF 00207 ! B = 1 - c.c 00208 ZB = 1. - PCOEF*PCOEF 00209 DO JLOOP1=1, KN 00210 ZSINLA = SIN(PLAC(JLOOP1)) 00211 ZCOSLA = COS(PLAC(JLOOP1)) 00212 ZSINLO = SIN(PLOC(JLOOP1)) 00213 ZCOSLO = COS(PLOC(JLOOP1)) 00214 ! X = cos(Lp-L) 00215 ZX = ZCOSLO*ZCOSSTRETCHLO + ZSINLO*ZSINSTRETCHLO 00216 ! Y = sin(Lp-L) 00217 ZY = ZSINSTRETCHLO*ZCOSLO - ZSINLO*ZCOSSTRETCHLO 00218 ! D = (1+c.c) + (1-c.c)(sin lp.sin l + cos lp.cos l.cos(Lp-L)) 00219 ZD = ZA + ZB*(ZSINSTRETCHLA*ZSINLA+ZCOSSTRETCHLA*ZCOSLA*ZX) 00220 ! (1-c.c)+(1+c.c)((sin lp.sin l + cos lp.cos l.cos(Lp-L)) 00221 ! sin lr = ------------------------------------------------------- 00222 ! D 00223 ZSINLAS = (ZB + ZA*(ZSINSTRETCHLA*ZSINLA+ZCOSSTRETCHLA*ZCOSLA*ZX)) / ZD 00224 ! D' = D * cos lr 00225 ZD = ZD * (AMAX1(1e-6,SQRT(1.-ZSINLAS*ZSINLAS))) 00226 ! 2.c.(cos lp.sin l - sin lp.cos l.cos(Lp-L)) 00227 ! cos Lr = ------------------------------------------- 00228 ! D' 00229 ZCOSLOS = 2.*PCOEF*(ZCOSSTRETCHLA*ZSINLA-ZSINSTRETCHLA*ZCOSLA*ZX) / ZD 00230 ! 2.c.cos l.cos(Lp-L) 00231 ! sin Lr = ------------------- 00232 ! D' 00233 ZSINLOS = 2.*PCOEF*(ZCOSLA*ZY) / ZD 00234 ! saturations (corrects calculation errors) 00235 ZSINLAS = MAX(ZSINLAS,-1.) 00236 ZSINLAS = MIN(ZSINLAS, 1.) 00237 ZCOSLOS = MAX(ZCOSLOS,-1.) 00238 ZCOSLOS = MIN(ZCOSLOS, 1.) 00239 ! back from sine & cosine 00240 PLAC(JLOOP1) = ASIN(ZSINLAS) 00241 IF (ZSINLOS>0) THEN 00242 PLOC(JLOOP1) = ACOS(ZCOSLOS) 00243 ELSE 00244 PLOC(JLOOP1) = -ACOS(ZCOSLOS) 00245 ENDIF 00246 ENDDO 00247 PLOC(:) = PLOC(:) * 180. / XPI 00248 PLAC(:) = PLAC(:) * 180. / XPI 00249 IF (LHOOK) CALL DR_HOOK('ARPEGE_STRETCH_A',1,ZHOOK_HANDLE) 00250 END SUBROUTINE ARPEGE_STRETCH_A 00251 !------------------------------------------------------------------------------- 00252 END SUBROUTINE HOR_INTERPOL_GAUSS
1.8.0
