SURFEX v7.3
General documentation of Surfex
|
00001 !------------------------------------------------------------------------------- 00002 ! ############################################################# 00003 SUBROUTINE PGD_TOPD(HPROGRAM) 00004 ! ############################################################# 00005 ! 00006 !!**** *PGD_TOPD* - routine to determine the masks that permit to couple ISBA grid with Topmodel one 00007 !! 00008 !! PURPOSE 00009 !! ------- 00010 !! 00011 !!** METHOD 00012 !! ------ 00013 !! 00014 !! EXTERNAL 00015 !! -------- 00016 !! 00017 !! IMPLICIT ARGUMENTS 00018 !! ------------------ 00019 !! 00020 !! REFERENCE 00021 !! --------- 00022 !! (from initial version of init_coupl_topo.f90 by K. Chancibault, modified by 00023 !! M. Lelay and B. Vincendon) 00024 !! 00025 !! AUTHOR 00026 !! ------ 00027 !! B. Vincendon *Meteo France* 00028 !! 00029 !! MODIFICATIONS 00030 !! ------------- 00031 !! Original 11/2011 00032 !------------------------------------------------------------------------------- 00033 ! 00034 !* 0. DECLARATIONS 00035 ! ------------ 00036 ! 00037 USE MODD_TOPODYN, ONLY : CCAT, NNCAT, XRTOP_D2, NMESHT, XDXT 00038 USE MODD_COUPLING_TOPD, ONLY : LCOUPL_TOPD, NIMAX, NJMAX, & 00039 XXI, XYI, NMASKI, NMASKT, NNPIX,& 00040 NNBV_IN_MESH, XBV_IN_MESH, XTOTBV_IN_MESH 00041 USE MODD_DUMMY_EXP_PROFILE, ONLY : XF_PARAM_BV, XC_DEPTH_RATIO_BV 00042 ! 00043 USE MODD_SURF_PAR, ONLY : NUNDEF 00044 USE MODD_SURF_ATM_GRID_N, ONLY : XGRID_PAR, CGRID !mll: ajout CGRID 00045 USE MODD_SURF_ATM_n, ONLY : NDIM_FULL 00046 ! 00047 USE MODD_ISBA_n, ONLY : CISBA 00048 ! 00049 USE MODE_GRIDTYPE_CONF_PROJ 00050 USE MODE_GRIDTYPE_LONLAT_REG 00051 USE MODE_GRIDTYPE_IGN 00052 ! 00053 USE MODI_GET_LUOUT 00054 USE MODI_READ_NAM_PGD_TOPD 00055 USE MODI_INIT_TOPD 00056 USE MODI_ABOR1_SFX 00057 USE MODI_MAKE_MASK_TOPD_TO_ISBA 00058 USE MODI_MAKE_MASK_ISBA_TO_TOPD 00059 USE MODI_WRITE_FILE_MASKTOPD 00060 USE MODI_OPEN_FILE 00061 USE MODI_CLOSE_FILE 00062 USE MODI_TOPD_TO_ISBA_SLOPE 00063 ! 00064 USE YOMHOOK ,ONLY : LHOOK, DR_HOOK 00065 USE PARKIND1 ,ONLY : JPRB 00066 ! 00067 IMPLICIT NONE 00068 ! 00069 !* 0.1 Declarations of arguments 00070 ! ------------------------- 00071 ! 00072 CHARACTER(LEN=*), INTENT(IN) :: HPROGRAM ! 00073 ! 00074 CHARACTER(LEN=50),DIMENSION(NNCAT) :: CNAME 00075 INTEGER :: IL ! number of points 00076 INTEGER :: JJ,JI,JK,JWRK ! loop control 00077 INTEGER :: JCAT,JMESH,JPIX ! loop control 00078 INTEGER :: ILUOUT ! Logical unit for output filr 00079 INTEGER :: IUNIT ! Logical unit for carte_asat file 00080 INTEGER :: IMESHL ! number of ISBA grid nodes 00081 ! 00082 REAL, DIMENSION(:), ALLOCATABLE :: ZXI, ZYI ! natural coordinates of ISBA grid (conformal projection or latlon) 00083 REAL, DIMENSION(:), ALLOCATABLE :: ZDXI, ZDYI ! Isba grid resolution in the conformal projection 00084 REAL, DIMENSION(:), ALLOCATABLE :: ZXN, ZYN ! isba nodes coordinates in the Lambert II coordinates 00085 REAL, DIMENSION(:), ALLOCATABLE :: ZLAT,ZLON ! Isba nodes geographical coordinates 00086 REAL, DIMENSION(:), ALLOCATABLE :: ZDTAV ! Averaged depth soil on TOP-LAT grid 00087 REAL :: ZLAT0 ! reference latitude 00088 REAL :: ZLON0 ! reference longitude 00089 REAL :: ZLONMIN,ZLONMAX ! min and max longitude values (latlon coordinates) 00090 REAL :: ZLATMIN,ZLATMAX ! min and max latitude values (latlon coordinates) 00091 REAL :: ZRPK ! projection parameter 00092 ! ! K=1 : stereographic north pole 00093 ! ! 0<K<1 : Lambert, north hemisphere 00094 ! ! K=0 : Mercator 00095 ! !-1<K<0 : Lambert, south hemisphere 00096 ! ! K=-1: stereographic south pole 00097 REAL :: ZBETA ! angle between grid and reference longitude 00098 REAL :: ZLATOR ! latitude of point of coordinates X=0, Y=0 00099 REAL :: ZLONOR ! longitude of point of coordinates X=0, Y=0! 00100 REAL, DIMENSION(:), ALLOCATABLE :: ZF_PARAM,ZC_DEPTH_RATIO ! 00101 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00102 !------------------------------------------------------------------------------- 00103 IF (LHOOK) CALL DR_HOOK('PGD_TOPD',0,ZHOOK_HANDLE) 00104 ! 00105 CALL GET_LUOUT(HPROGRAM,ILUOUT) 00106 ! 00107 CALL READ_NAM_PGD_TOPD(HPROGRAM,LCOUPL_TOPD,CCAT,XF_PARAM_BV,XC_DEPTH_RATIO_BV) 00108 ! 00109 IF (LCOUPL_TOPD .AND. CISBA/='3-L') & 00110 CALL ABOR1_SFX("PGD_TOPD: coupling with topmodel only runs with CISBA=3-L") 00111 ! 00112 ! 1. Reads the namelists 00113 ! -------------------- 00114 IF (LCOUPL_TOPD) THEN 00115 ! 00116 WRITE(ILUOUT,*) 'Debut pgd_topd' 00117 ! 00118 ! 3. Initialises variables specific to Topmodel 00119 ! ------------------------------------------- 00120 WRITE(ILUOUT,*) 'NNCAT',NNCAT 00121 ! 00122 CALL INIT_TOPD(HPROGRAM) 00123 ! 00124 ! 4. Compute masks to couple ISBA and TOPODYN grids 00125 ! ------------------------------------------------------- 00126 ! 00127 !* 1. Calcul of the plane coordinates 00128 ! ------------------------------- 00129 ! 00130 !* 1.1 Gestion des projections conformes 00131 ! --------------------------------- 00132 ! 00133 ALLOCATE(NMASKT(NNCAT,NMESHT)) 00134 NMASKT(:,:)=NUNDEF 00135 ! 00136 IF(CGRID.EQ.'CONF PROJ') THEN 00137 ! 00138 WRITE(ILUOUT,*) 'GRILLE PROJ CONF (application Cevennes)' 00139 ! 00140 !* 1.1.1 lecture des coordonnees X et Y conformes 00141 ! ----------------------------------------- 00142 ! 00143 ALLOCATE(ZXI(NDIM_FULL)) 00144 ALLOCATE(ZYI(NDIM_FULL)) 00145 ZXI(:)=0.0 00146 ZYI(:)=0.0 00147 ! 00148 ALLOCATE(ZDXI(NDIM_FULL)) 00149 ALLOCATE(ZDYI(NDIM_FULL)) 00150 ! 00151 CALL GET_GRIDTYPE_CONF_PROJ(XGRID_PAR,PLAT0=ZLAT0,PLON0=ZLON0,PRPK=ZRPK, & 00152 PBETA=ZBETA,PLATOR=ZLATOR,PLONOR=ZLONOR, & 00153 KIMAX=NIMAX,KJMAX=NJMAX,PX=ZXI,PY=ZYI, & 00154 PDX=ZDXI,PDY=ZDYI) 00155 ! 00156 IMESHL = (NIMAX+1)*(NJMAX+1) 00157 ! 00158 !* 1.1.2 calcul des coordonnees X et Y conformes des noeuds de la grille ISBA 00159 ! -------------------------------------------------------------------- 00160 ALLOCATE(ZXN(IMESHL)) 00161 ALLOCATE(ZYN(IMESHL)) 00162 ! 00163 DO JJ=1,NIMAX 00164 ZXN(JJ) = ZXI(JJ) - ZDXI(JJ)/2. 00165 ENDDO 00166 ZXN(NIMAX+1) = ZXI(NIMAX) + ZDXI(NIMAX)/2. 00167 ! 00168 DO JJ=1,NJMAX 00169 JWRK = (JJ-1)*(NIMAX+1)+1 ! indice sur la grille des noeuds 00170 JI = (JJ-1)*NIMAX+1 ! indice sur la grille des mailles 00171 ZYN(JWRK) = ZYI(JI) - ZDYI(JI)/2. 00172 ENDDO 00173 ! 00174 JJ = ((NJMAX+1)-1)*(NIMAX+1)+1 ! Indice sur la grille des noeuds 00175 JI = (NJMAX-1)*NIMAX+1 00176 ZYN(JJ) = ZYI(JI) + ZDYI(JI)/2. 00177 ! 00178 DEALLOCATE(ZDXI) 00179 DEALLOCATE(ZDYI) 00180 DEALLOCATE(ZXI) 00181 DEALLOCATE(ZYI) 00182 ! 00183 DO JJ=1,NIMAX+1 00184 DO JI=2,NJMAX+1 00185 JK = (JI-1)*(NIMAX+1)+JJ 00186 ZXN(JK) = ZXN(JJ) 00187 ENDDO 00188 ENDDO 00189 ! 00190 DO JI=1,NJMAX+1 00191 DO JJ=2,NIMAX+1 00192 JK = (JI-1)*(NIMAX+1)+JJ 00193 JWRK = (JI-1)*(NIMAX+1)+1 00194 ZYN(JK) = ZYN(JWRK) 00195 ENDDO 00196 ENDDO 00197 ! 00198 !* 1.1.3 calcul des coordonnées géographiques des noeuds de la grille ISBA 00199 ! ----------------------------------------------------------------- 00200 ALLOCATE(ZLAT(IMESHL)) 00201 ALLOCATE(ZLON(IMESHL)) 00202 CALL LATLON_CONF_PROJ(ZLAT0,ZLON0,ZRPK,ZBETA,ZLATOR,ZLONOR,ZXN,ZYN,ZLAT,ZLON) 00203 DEALLOCATE(ZXN) 00204 DEALLOCATE(ZYN) 00205 ! 00206 !* 1.2 Gestion des coordonnees geographiques 00207 ! ------------------------------------- 00208 ! 00209 ELSE IF(CGRID.EQ.'LONLAT REG') THEN 00210 ! 00211 WRITE(ILUOUT,*) 'GRILLE LONLAT REG (application AMMA)' 00212 ! 00213 ALLOCATE(ZXI(NDIM_FULL)) 00214 ALLOCATE(ZYI(NDIM_FULL)) 00215 ZXI(:)=0.0 00216 ZYI(:)=0.0 00217 CALL GET_GRIDTYPE_LONLAT_REG(XGRID_PAR,PLONMIN=ZLONMIN,PLONMAX=ZLONMAX, & 00218 PLATMIN=ZLATMIN,PLATMAX=ZLATMAX,KLON=NIMAX,KLAT=NJMAX, & 00219 KL=IL,PLON=ZXI,PLAT=ZYI) 00220 ! 00221 IMESHL=(NIMAX+1)*(NJMAX+1) 00222 ! 00223 ALLOCATE(ZLON(IMESHL)) 00224 ALLOCATE(ZLAT(IMESHL)) 00225 ALLOCATE(ZDXI(NDIM_FULL)) 00226 ALLOCATE(ZDYI(NDIM_FULL)) 00227 ! 00228 ZDXI(:)=(ZLONMAX-ZLONMIN)/(NIMAX-1) 00229 ZDYI(:)=(ZLATMAX-ZLATMIN)/(NJMAX-1) 00230 ! 00231 DO JJ=1,NIMAX 00232 ZLON(JJ) = ZXI(JJ) - ZDXI(JJ)/2. 00233 ENDDO 00234 ZLON(NIMAX+1) = ZXI(NIMAX) + ZDXI(NIMAX)/2. 00235 ! 00236 DO JJ=1,NJMAX 00237 JWRK=(JJ-1)*(NIMAX+1)+1 ! indice sur la grille des noeuds 00238 JI=(JJ-1)*NIMAX+1 ! indice sur la grille des mailles 00239 ZLAT(JWRK) = ZYI(JI) - ZDYI(JI)/2. 00240 ENDDO 00241 ! 00242 JJ=((NJMAX+1)-1)*(NIMAX+1)+1 ! Indice sur la grille des noeuds 00243 JI=(NJMAX-1)*NIMAX+1 00244 ZLAT(JJ) = ZYI(JI) + ZDYI(JI)/2. 00245 ! 00246 DEALLOCATE(ZDXI) 00247 DEALLOCATE(ZDYI) 00248 DEALLOCATE(ZXI) 00249 DEALLOCATE(ZYI) 00250 ! 00251 DO JJ=1,NIMAX+1 00252 DO JI=2,NJMAX+1 00253 JK = (JI-1)*(NIMAX+1)+JJ 00254 ZLON(JK) = ZLON(JJ) 00255 ENDDO 00256 ENDDO 00257 ! 00258 DO JI=1,NJMAX+1 00259 DO JJ=2,NIMAX+1 00260 JK=(JI-1)*(NIMAX+1)+JJ 00261 JWRK=(JI-1)*(NIMAX+1)+1 00262 ZLAT(JK)=ZLAT(JWRK) 00263 ENDDO 00264 ENDDO 00265 ! 00266 ELSE 00267 ! 00268 WRITE(ILUOUT,*) 'ERREUR: TYPE DE GRILLE NON GERE PAR LE CODE' 00269 CALL ABOR1_SFX("PGD_TOPD: TYPE DE GRILLE NON GERE PAR LE CODE") 00270 ! 00271 ENDIF 00272 ! 00273 !* 2.0 calcul des coordonnées lambert II étendu des noeuds de la grille ISBA 00274 ! ---------------------------------------------------------------------- 00275 ! 00276 ALLOCATE(XXI(IMESHL)) 00277 ALLOCATE(XYI(IMESHL)) 00278 CALL XY_IGN(5,XXI,XYI,ZLAT,ZLON) 00279 DEALLOCATE(ZLAT) 00280 DEALLOCATE(ZLON) 00281 ! 00282 !* 2.0 mask 00283 ! ---- 00284 !*** 00285 CALL MAKE_MASK_TOPD_TO_ISBA(NDIM_FULL) 00286 ! 00287 ALLOCATE(NNPIX(NDIM_FULL)) 00288 NNPIX(:) = NUNDEF 00289 DO JJ=1,NDIM_FULL 00290 NNPIX(JJ) = COUNT(NMASKT(:,:)==JJ) 00291 ENDDO 00292 ! 00293 CALL MAKE_MASK_ISBA_TO_TOPD(NDIM_FULL) 00294 ! 00295 CALL WRITE_FILE_MASKTOPD(NDIM_FULL) 00296 ! 00297 !* 3.0 Compute Mean slope over each ISBA_MESH 00298 ! ---------------------------------------------------------------------- 00299 CALL TOPD_TO_ISBA_SLOPE(NDIM_FULL) 00300 ! 00301 !* 4.0 Compute F and DC for each ISBA mesh 00302 ! ---------------------------------------------------------------------- 00303 ! 00304 ALLOCATE(NNBV_IN_MESH(NDIM_FULL,NNCAT)) 00305 ALLOCATE(XBV_IN_MESH(NDIM_FULL,NNCAT)) 00306 ALLOCATE(XTOTBV_IN_MESH(NDIM_FULL)) 00307 ! 00308 XTOTBV_IN_MESH(:) = 0.0 00309 ! 00310 DO JMESH=1,NDIM_FULL 00311 XBV_IN_MESH(JMESH,:)=0.0 00312 DO JCAT=1,NNCAT 00313 NNBV_IN_MESH(JMESH,JCAT) = COUNT(NMASKI(JMESH,JCAT,:)/=NUNDEF) 00314 XBV_IN_MESH(JMESH,JCAT) = REAL(NNBV_IN_MESH(JMESH,JCAT))*XDXT(JCAT)**2 00315 XTOTBV_IN_MESH(JMESH) = XTOTBV_IN_MESH(JMESH) + XBV_IN_MESH(JMESH,JCAT) 00316 ENDDO 00317 ENDDO 00318 ! 00319 ALLOCATE (ZF_PARAM(NDIM_FULL)) 00320 ALLOCATE (ZC_DEPTH_RATIO(NDIM_FULL)) 00321 ! 00322 ZF_PARAM(:) = 0. 00323 ZC_DEPTH_RATIO(:) = 0. 00324 DO JCAT=1,NNCAT 00325 DO JMESH=1,NDIM_FULL 00326 IF ( XTOTBV_IN_MESH(JMESH)/=0. ) THEN 00327 ZF_PARAM(JMESH) = ZF_PARAM(JMESH) + XF_PARAM_BV(JCAT)*XBV_IN_MESH(JMESH,JCAT)/XTOTBV_IN_MESH(JMESH) 00328 ZC_DEPTH_RATIO(JMESH) = ZC_DEPTH_RATIO(JMESH) + XC_DEPTH_RATIO_BV(JCAT)*XBV_IN_MESH(JMESH,JCAT)/XTOTBV_IN_MESH(JMESH) 00329 ENDIF 00330 ENDDO 00331 ENDDO 00332 ! 00333 WHERE (ZF_PARAM==0.) 00334 ZF_PARAM=2.5 00335 ZC_DEPTH_RATIO=1. 00336 ENDWHERE 00337 ! 00338 !write(*,*) 'f min max isba',MINVAL(ZF_PARAM),MAXVAL(ZF_PARAM) 00339 !write(*,*) 'dc min max isba',MINVAL(ZC_DEPTH_RATIO),MAXVAL(ZC_DEPTH_RATIO) 00340 ! 00341 CALL OPEN_FILE('ASCII ',IUNIT,'carte_f_dc.txt','FORMATTED',HACTION='WRITE') 00342 DO JMESH=1,NDIM_FULL 00343 WRITE(IUNIT,*) ZF_PARAM(JMESH),ZC_DEPTH_RATIO(JMESH) 00344 ENDDO 00345 CALL CLOSE_FILE('ASCII ',IUNIT) 00346 ! 00347 DEALLOCATE(ZF_PARAM) 00348 DEALLOCATE(ZC_DEPTH_RATIO) 00349 ! 00350 WRITE(ILUOUT,*) 'Couplage avec TOPMODEL active' 00351 ! 00352 ELSE 00353 ! 00354 WRITE(ILUOUT,*) 'Pas de couplage avec TOPMODEL' 00355 ! 00356 ENDIF 00357 ! 00358 IF (LHOOK) CALL DR_HOOK('PGD_TOPD',1,ZHOOK_HANDLE) 00359 ! 00360 END SUBROUTINE PGD_TOPD