|
SURFEX v7.3
General documentation of Surfex
|
|
SURFEX v7.3
General documentation of Surfex
|
00001 ! ##################### 00002 MODULE MODE_READ_GRIB 00003 ! ##################### 00004 !------------------------------------------------------------------- 00005 ! 00006 USE MODI_ABOR1_SFX 00007 USE GRIB_API 00008 ! 00009 USE YOMHOOK ,ONLY : LHOOK, DR_HOOK 00010 USE PARKIND1 ,ONLY : JPRB 00011 ! 00012 CONTAINS 00013 !------------------------------------------------------------------- 00014 ! #################### 00015 SUBROUTINE MAKE_GRIB_INDEX(HGRIB) 00016 ! #################### 00017 ! 00018 USE MODD_GRID_GRIB, ONLY : CGRIB_FILE, NIDX 00019 ! 00020 IMPLICIT NONE 00021 ! 00022 CHARACTER(LEN=*), INTENT(IN) :: HGRIB 00023 ! 00024 INTEGER(KIND=kindOfInt) :: IRET 00025 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00026 ! 00027 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:MAKE_GRIB_INDEX',0,ZHOOK_HANDLE) 00028 ! 00029 IF (CGRIB_FILE==HGRIB) CALL DR_HOOK('MODE_READ_GRIB:MAKE_GRIB_INDEX',1,ZHOOK_HANDLE) 00030 IF (CGRIB_FILE==HGRIB) RETURN 00031 ! 00032 CGRIB_FILE=HGRIB 00033 ! 00034 CALL GRIB_INDEX_CREATE(NIDX,HGRIB,'indicatorOfParameter',IRET) 00035 IF (IRET/=0) CALL ABOR1_SFX("MODE_READ_GRIB:MAKE_GRIB_INDEX: error while creating the grib index") 00036 ! 00037 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:MAKE_GRIB_INDEX',1,ZHOOK_HANDLE) 00038 ! 00039 END SUBROUTINE MAKE_GRIB_INDEX 00040 !------------------------------------------------------------------- 00041 ! #################### 00042 SUBROUTINE CLEAR_GRIB_INDEX 00043 ! #################### 00044 ! 00045 USE MODD_GRID_GRIB, ONLY : CGRIB_FILE, NIDX 00046 ! 00047 IMPLICIT NONE 00048 ! 00049 INTEGER(KIND=kindOfInt) :: IRET 00050 ! 00051 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00052 ! 00053 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:CLEAR_GRIB_INDEX',0,ZHOOK_HANDLE) 00054 ! 00055 IF (CGRIB_FILE.NE."") THEN 00056 CGRIB_FILE="" 00057 CALL GRIB_INDEX_RELEASE(NIDX,IRET) 00058 IF (IRET/=0) CALL ABOR1_SFX("MODE_READ_GRIB:MAKE_GRIB_INDEX: error while deleting the grib index") 00059 ENDIF 00060 ! 00061 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:CLEAR_GRIB_INDEX',1,ZHOOK_HANDLE) 00062 ! 00063 END SUBROUTINE CLEAR_GRIB_INDEX 00064 !------------------------------------------------------------------- 00065 ! #################### 00066 SUBROUTINE GET_GRIB_MESSAGE(KLUOUT,KLTYPE,KLEV1,KLEV2,KGRIB,KFOUND) 00067 ! #################### 00068 ! 00069 USE MODD_GRID_GRIB, ONLY : NIDX 00070 ! 00071 IMPLICIT NONE 00072 ! 00073 INTEGER, INTENT(IN) :: KLUOUT 00074 INTEGER, INTENT(INOUT) :: KLTYPE 00075 INTEGER, INTENT(INOUT) :: KLEV1 00076 INTEGER, INTENT(INOUT) :: KLEV2 00077 INTEGER(KIND=kindOfInt), INTENT(INOUT) :: KGRIB 00078 INTEGER, INTENT(OUT) :: KFOUND 00079 ! 00080 INTEGER :: ILTYPE 00081 INTEGER :: ILEV1 00082 INTEGER :: ILEV2 00083 INTEGER(KIND=kindOfInt) :: IRET 00084 ! 00085 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00086 ! 00087 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:GET_GRIB_MESSAGE',0,ZHOOK_HANDLE) 00088 ! 00089 IRET = 0 00090 KFOUND=0 00091 ! 00092 DO WHILE (IRET /= GRIB_END_OF_INDEX .AND. KFOUND/=3) 00093 ! 00094 IRET = 0 00095 KFOUND=0 00096 ! 00097 IF (KLTYPE/=-2) THEN 00098 CALL GRIB_GET(KGRIB,'indicatorOfTypeOfLevel',ILTYPE,IRET) 00099 CALL TEST_IRET(KLUOUT,ILTYPE,KLTYPE,IRET) 00100 ENDIF 00101 ! 00102 IF (IRET.EQ.0) THEN 00103 ! 00104 KFOUND = KFOUND + 1 00105 ! 00106 IF (KLEV1/=-2) THEN 00107 CALL GRIB_GET(KGRIB,'topLevel',ILEV1,IRET) 00108 CALL TEST_IRET(KLUOUT,ILEV1,KLEV1,IRET) 00109 ENDIF 00110 ! 00111 IF (IRET.EQ.0) THEN 00112 ! 00113 KFOUND = KFOUND + 1 00114 ! 00115 IF (KLEV2/=-2) THEN 00116 CALL GRIB_GET(KGRIB,'bottomLevel',ILEV2,IRET) 00117 CALL TEST_IRET(KLUOUT,ILEV2,KLEV2,IRET) 00118 ENDIF 00119 ! 00120 IF (IRET.EQ.0) KFOUND = KFOUND + 1 00121 ! 00122 ENDIF 00123 ! 00124 ENDIF 00125 ! 00126 IF (KFOUND.NE.3) THEN 00127 CALL GRIB_RELEASE(KGRIB) 00128 CALL GRIB_NEW_FROM_INDEX(NIDX,KGRIB,IRET) 00129 ENDIF 00130 ! 00131 ENDDO 00132 ! 00133 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:GET_GRIB_MESSAGE',1,ZHOOK_HANDLE) 00134 ! 00135 CONTAINS 00136 ! 00137 ! ############## 00138 SUBROUTINE TEST_IRET(KLUOUT,VAL1,VAL0,KRET) 00139 ! ############## 00140 ! 00141 IMPLICIT NONE 00142 ! 00143 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00144 INTEGER, INTENT(IN) :: VAL1 00145 INTEGER, INTENT(INOUT) :: VAL0 00146 INTEGER(KIND=kindOfInt), INTENT(INOUT) :: KRET ! number of the message researched 00147 ! 00148 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00149 ! 00150 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:TEST_IRET',0,ZHOOK_HANDLE) 00151 ! 00152 IF (KRET > 0) THEN 00153 WRITE (KLUOUT,'(A)')' | Error encountered in the Grib file, skipping field' 00154 ELSE IF (KRET == -6) THEN 00155 WRITE (KLUOUT,'(A)')' | ECMWF pseudo-Grib data encountered, skipping field' 00156 ELSEIF (VAL1 /= VAL0) THEN 00157 IF (VAL0 == -1) THEN 00158 VAL0 = VAL1 00159 ELSE 00160 KRET=1 00161 ENDIF 00162 ENDIF 00163 ! 00164 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:TEST_IRET',1,ZHOOK_HANDLE) 00165 END SUBROUTINE TEST_IRET 00166 ! 00167 END SUBROUTINE GET_GRIB_MESSAGE 00168 !------------------------------------------------------------------- 00169 ! #################### 00170 SUBROUTINE READ_GRIB(HGRIB,KLUOUT,KPARAM,KRET,PFIELD,KLTYPE,KLEV1,KLEV2,KPARAM2) 00171 ! #################### 00172 ! 00173 USE MODD_GRID_GRIB, ONLY : NIDX 00174 ! 00175 IMPLICIT NONE 00176 ! 00177 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! name of the GRIB file 00178 INTEGER, INTENT(IN) :: KLUOUT 00179 INTEGER,INTENT(IN) :: KPARAM ! Parameter to read 00180 INTEGER(KIND=kindOfInt), INTENT(OUT) :: KRET 00181 REAL, DIMENSION(:), POINTER :: PFIELD 00182 INTEGER,INTENT(INOUT), OPTIONAL :: KLTYPE ! Level type 00183 INTEGER,INTENT(INOUT), OPTIONAL :: KLEV1 ! Level parameter 1 00184 INTEGER,INTENT(INOUT), OPTIONAL :: KLEV2 ! Level parameter 2 00185 INTEGER, INTENT(INOUT), OPTIONAL :: KPARAM2 00186 ! 00187 INTEGER :: ILTYPE, ILEV1, ILEV2 00188 INTEGER(KIND=kindOfInt) :: IGRIB 00189 INTEGER :: ISIZE, IFOUND 00190 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00191 ! 00192 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB',0,ZHOOK_HANDLE) 00193 ! 00194 ILTYPE=-2 00195 IF (PRESENT(KLTYPE)) ILTYPE=KLTYPE 00196 ILEV1=-2 00197 IF (PRESENT(KLEV1)) ILEV1=KLEV1 00198 ILEV2=-2 00199 IF (PRESENT(KLEV2)) ILEV2=KLEV2 00200 ! 00201 CALL MAKE_GRIB_INDEX(HGRIB) 00202 ! 00203 IFOUND=0 00204 KRET=0 00205 ! 00206 CALL GRIB_INDEX_SELECT(NIDX,'indicatorOfParameter',KPARAM,KRET) 00207 CALL GRIB_NEW_FROM_INDEX(NIDX,IGRIB,KRET) 00208 IF (KRET.EQ.0) CALL GET_GRIB_MESSAGE(KLUOUT,ILTYPE,ILEV1,ILEV2,IGRIB,IFOUND) 00209 ! 00210 IF (PRESENT(KPARAM2)) THEN 00211 IF (IFOUND/=3) THEN 00212 CALL GRIB_INDEX_SELECT(NIDX,'indicatorOfParameter',KPARAM2,KRET) 00213 CALL GRIB_NEW_FROM_INDEX(NIDX,IGRIB,KRET) 00214 IF (KRET.EQ.0) THEN 00215 ILTYPE=-2 00216 IF (PRESENT(KLTYPE)) ILTYPE=KLTYPE 00217 CALL GET_GRIB_MESSAGE(KLUOUT,ILTYPE,ILEV1,ILEV2,IGRIB,IFOUND) 00218 ENDIF 00219 ELSE 00220 KPARAM2 = KPARAM 00221 ENDIF 00222 ENDIF 00223 ! 00224 IF (IFOUND==3) THEN 00225 ! 00226 IF (PRESENT(KLTYPE)) KLTYPE = ILTYPE 00227 IF (PRESENT(KLEV1)) KLEV1 = ILEV1 00228 IF (PRESENT(KLEV2)) KLEV2 = ILEV2 00229 ! 00230 IF (.NOT.ASSOCIATED(PFIELD)) THEN 00231 CALL GRIB_GET_SIZE(IGRIB,'values',ISIZE,KRET) 00232 IF (KRET.NE.0) CALL ABOR1_SFX("MODE_READ_GRIB:READ_GRIB: Problem getting size of values") 00233 ALLOCATE(PFIELD(ISIZE)) 00234 ENDIF 00235 ! 00236 CALL GRIB_GET(IGRIB,'values',PFIELD,KRET) 00237 IF (KRET.NE.0) CALL ABOR1_SFX("MODE_READ_GRIB:READ_GRIB: Problem getting values") 00238 CALL GRIB_RELEASE(IGRIB,KRET) 00239 IF (KRET.NE.0) CALL ABOR1_SFX("MODE_READ_GRIB:READ_GRIB: Problem releasing memory") 00240 ! 00241 ELSE 00242 ! 00243 KRET = 1 00244 ! 00245 ENDIF 00246 ! 00247 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB',1,ZHOOK_HANDLE) 00248 END SUBROUTINE READ_GRIB 00249 ! 00250 !------------------------------------------------------------------- 00251 !------------------------------------------------------------------- 00252 ! #################### 00253 SUBROUTINE READ_GRIB_LAND_MASK(HGRIB,KLUOUT,HINMODEL,PMASK) 00254 ! #################### 00255 ! 00256 IMPLICIT NONE 00257 ! 00258 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00259 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00260 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00261 REAL, DIMENSION(:), POINTER :: PMASK ! Land mask 00262 ! 00263 INTEGER(KIND=kindOfInt) :: IRET ! return code 00264 INTEGER :: ILTYPE ! leveltype 00265 INTEGER :: ILEV ! level 00266 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00267 !------------------------------------------------------------------- 00268 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_LAND_MASK',0,ZHOOK_HANDLE) 00269 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_LAND_MASK: | Reading land mask from ',HINMODEL 00270 ! 00271 SELECT CASE (HINMODEL) 00272 CASE ('ECMWF ') 00273 CALL READ_GRIB(HGRIB,KLUOUT,172,IRET,PMASK) 00274 CASE ('ARPEGE','ALADIN','MOCAGE') 00275 CALL READ_GRIB(HGRIB,KLUOUT,81,IRET,PMASK) 00276 CASE ('HIRLAM') 00277 ILTYPE=105 00278 ILEV =0 00279 CALL READ_GRIB(HGRIB,KLUOUT,81,IRET,PMASK,KLTYPE=ILTYPE,KLEV1=ILEV) 00280 CASE DEFAULT 00281 CALL ABOR1_SFX('MODE_READ_GRIB:READ_GRIB_LAND_MASK: OPTION NOT SUPPORTED '//HINMODEL) 00282 END SELECT 00283 ! 00284 IF (IRET /= 0) THEN 00285 CALL ABOR1_SFX('MODE_READ_GRIB: LAND SEA MASK MISSING (READ_GRIB_LAND_MASK)') 00286 END IF 00287 ! 00288 WHERE (PMASK>0.5) 00289 PMASK = 1. 00290 ELSEWHERE 00291 PMASK = 0. 00292 END WHERE 00293 ! 00294 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_LAND_MASK',1,ZHOOK_HANDLE) 00295 END SUBROUTINE READ_GRIB_LAND_MASK 00296 !------------------------------------------------------------------- 00297 ! ############################ 00298 SUBROUTINE READ_GRIB_ZS(HGRIB,KLUOUT,HINMODEL,PZS) 00299 ! ############################ 00300 ! 00301 USE MODD_CSTS, ONLY : XG 00302 ! 00303 IMPLICIT NONE 00304 ! 00305 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00306 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00307 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00308 REAL, DIMENSION(:), POINTER :: PZS ! 00309 ! 00310 INTEGER(KIND=kindOfInt) :: IRET ! return code 00311 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00312 !------------------------------------------------------------------- 00313 !* Read orography 00314 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_ZS',0,ZHOOK_HANDLE) 00315 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_ZS: | Reading orography from ',HINMODEL 00316 ! 00317 SELECT CASE (HINMODEL) 00318 CASE ('ECMWF ') 00319 CALL READ_GRIB(HGRIB,KLUOUT,129,IRET,PZS) 00320 CASE ('ARPEGE','MOCAGE') 00321 CALL READ_GRIB(HGRIB,KLUOUT,8,IRET,PZS) 00322 CASE ('HIRLAM','ALADIN') 00323 CALL READ_GRIB(HGRIB,KLUOUT,6,IRET,PZS) 00324 CASE DEFAULT 00325 CALL ABOR1_SFX('MODE_READ_GRIB:READ_GRIB_ZS:OPTION NOT SUPPORTED '//HINMODEL) 00326 END SELECT 00327 ! 00328 IF (IRET /= 0) THEN 00329 CALL ABOR1_SFX('MODE_READ_GRIB: OROGRAPHY MISSING (READ_GRIB_ZS_LAND)') 00330 END IF 00331 ! 00332 ! Datas given in archives are multiplied by the gravity acceleration 00333 PZS = PZS / XG 00334 ! 00335 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_ZS',1,ZHOOK_HANDLE) 00336 END SUBROUTINE READ_GRIB_ZS 00337 !------------------------------------------------------------------- 00338 ! ############################ 00339 SUBROUTINE READ_GRIB_ZS_LAND(HGRIB,KLUOUT,HINMODEL,PMASK,PZSL) 00340 ! ############################ 00341 ! 00342 IMPLICIT NONE 00343 ! 00344 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00345 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00346 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00347 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 00348 REAL, DIMENSION(:), POINTER :: PZSL ! 00349 ! 00350 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00351 !------------------------------------------------------------------- 00352 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_ZS_LAND',0,ZHOOK_HANDLE) 00353 ! 00354 CALL READ_GRIB_ZS(HGRIB,KLUOUT,HINMODEL,PZSL) 00355 ! 00356 IF (SIZE(PMASK)==SIZE(PZSL)) & 00357 WHERE (PMASK(:)/=1.) PZSL = 0. 00358 ! 00359 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_ZS_LAND',1,ZHOOK_HANDLE) 00360 END SUBROUTINE READ_GRIB_ZS_LAND 00361 !------------------------------------------------------------------- 00362 ! ############################ 00363 SUBROUTINE READ_GRIB_ZS_SEA(HGRIB,KLUOUT,HINMODEL,PMASK,PZSS) 00364 ! ############################ 00365 ! 00366 IMPLICIT NONE 00367 ! 00368 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00369 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00370 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00371 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 00372 REAL, DIMENSION(:), POINTER :: PZSS ! 00373 ! 00374 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00375 !------------------------------------------------------------------- 00376 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_ZS_SEA',0,ZHOOK_HANDLE) 00377 ! 00378 IF (HINMODEL=='HIRLAM') THEN 00379 CALL ABOR1_SFX('MODE_READ_GRIB:READ_GRIB_ZS_SEA:OPTION NOT SUPPORTED '//HINMODEL) 00380 ELSE 00381 CALL READ_GRIB_ZS(HGRIB,KLUOUT,HINMODEL,PZSS) 00382 ENDIF 00383 ! 00384 IF (SIZE(PMASK)==SIZE(PZSS)) & 00385 WHERE (PMASK(:)/=0.) PZSS = 0. 00386 ! 00387 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_ZS_SEA',1,ZHOOK_HANDLE) 00388 END SUBROUTINE READ_GRIB_ZS_SEA 00389 !------------------------------------------------------------------- 00390 ! ########################### 00391 SUBROUTINE READ_GRIB_T(HGRIB,KLUOUT,HINMODEL,PT) 00392 ! ########################### 00393 ! 00394 IMPLICIT NONE 00395 ! 00396 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00397 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00398 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00399 REAL, DIMENSION(:), POINTER :: PT ! 00400 ! 00401 INTEGER(KIND=kindOfInt) :: IRET ! return code 00402 INTEGER :: ILTYPE ! type of level (Grib code table 3) 00403 INTEGER :: ILEV ! level definition 00404 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00405 !------------------------------------------------------------------- 00406 !* Read surface temperature 00407 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_T',0,ZHOOK_HANDLE) 00408 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_T: | Reading surface temperature' 00409 ! 00410 SELECT CASE (HINMODEL) 00411 CASE ('ECMWF ') 00412 CALL READ_GRIB(HGRIB,KLUOUT,139,IRET,PT) 00413 00414 CASE ('ARPEGE','ALADIN','MOCAGE') 00415 ILEV=0 00416 ILTYPE=111 00417 CALL READ_GRIB(HGRIB,KLUOUT,11,IRET,PT,KLTYPE=ILTYPE,KLEV1=ILEV) 00418 IF (IRET /= 0) THEN 00419 ILTYPE=1 00420 CALL READ_GRIB(HGRIB,KLUOUT,11,IRET,PT,KLTYPE=ILTYPE) 00421 IF (IRET /= 0) THEN 00422 ILTYPE=105 00423 CALL READ_GRIB(HGRIB,KLUOUT,11,IRET,PT,KLTYPE=ILTYPE,KLEV1=ILEV) 00424 ENDIF 00425 END IF 00426 00427 CASE ('HIRLAM ') 00428 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_T: | Reading surface temperature tile 4' 00429 ILTYPE=105 00430 ILEV=904 00431 CALL READ_GRIB(HGRIB,KLUOUT,11,IRET,PT,KLTYPE=ILTYPE,KLEV1=ILEV) 00432 00433 CASE DEFAULT 00434 CALL ABOR1_SFX('MODE_READ_GRIB:READ_GRIB_T:OPTION NOT SUPPORTED '//HINMODEL) 00435 END SELECT 00436 ! 00437 IF (IRET /= 0) THEN 00438 CALL ABOR1_SFX('MODE_READ_GRIB: SURFACE TEMPERATURE MISSING (READ_GRIB_T)') 00439 END IF 00440 ! 00441 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_T',1,ZHOOK_HANDLE) 00442 END SUBROUTINE READ_GRIB_T 00443 !------------------------------------------------------------------- 00444 ! ########################### 00445 SUBROUTINE READ_GRIB_TS(HGRIB,KLUOUT,HINMODEL,PMASK,PTS) 00446 ! ########################### 00447 ! 00448 USE MODD_SURF_PAR, ONLY : XUNDEF 00449 ! 00450 IMPLICIT NONE 00451 ! 00452 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00453 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00454 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00455 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 00456 REAL, DIMENSION(:), POINTER :: PTS ! 00457 ! 00458 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00459 !------------------------------------------------------------------- 00460 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TS',0,ZHOOK_HANDLE) 00461 ! 00462 CALL READ_GRIB_T(HGRIB,KLUOUT,HINMODEL,PTS) 00463 ! 00464 IF (SIZE(PMASK)==SIZE(PTS)) & 00465 WHERE (PMASK(:)/=1.) PTS = XUNDEF 00466 ! 00467 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TS',1,ZHOOK_HANDLE) 00468 END SUBROUTINE READ_GRIB_TS 00469 !------------------------------------------------------------------- 00470 ! ########################### 00471 SUBROUTINE READ_GRIB_SST(HGRIB,KLUOUT,HINMODEL,PMASK,PSST) 00472 ! ########################### 00473 ! 00474 USE MODD_SURF_PAR, ONLY : XUNDEF 00475 ! 00476 IMPLICIT NONE 00477 ! 00478 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00479 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00480 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00481 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 00482 REAL, DIMENSION(:), POINTER :: PSST ! 00483 ! 00484 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00485 !------------------------------------------------------------------- 00486 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_SST',0,ZHOOK_HANDLE) 00487 ! 00488 IF (HINMODEL=='HIRLAM') THEN 00489 CALL ABOR1_SFX('MODE_READ_GRIB:READ_GRIB_SST:OPTION NOT SUPPORTED '//HINMODEL) 00490 ELSE 00491 CALL READ_GRIB_T(HGRIB,KLUOUT,HINMODEL,PSST) 00492 ENDIF 00493 ! 00494 IF (SIZE(PMASK)==SIZE(PSST)) & 00495 WHERE (PMASK(:)/=0.) PSST = XUNDEF 00496 ! 00497 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_SST',1,ZHOOK_HANDLE) 00498 END SUBROUTINE READ_GRIB_SST 00499 !------------------------------------------------------------------- 00500 ! ########################### 00501 SUBROUTINE READ_GRIB_T2(HGRIB,KLUOUT,HINMODEL,PMASK,PT2) 00502 ! ########################### 00503 ! 00504 USE MODD_SURF_PAR, ONLY : XUNDEF 00505 ! 00506 IMPLICIT NONE 00507 ! 00508 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00509 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00510 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00511 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 00512 REAL, DIMENSION(:), POINTER :: PT2 ! 00513 ! 00514 INTEGER(KIND=kindOfInt) :: IRET 00515 INTEGER :: ILTYPE ! type of level (Grib code table 3) 00516 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00517 !------------------------------------------------------------------- 00518 !* Read deep soil temperature 00519 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_T2',0,ZHOOK_HANDLE) 00520 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_T2: | Reading deep soil temperature' 00521 ! 00522 SELECT CASE (HINMODEL) 00523 CASE ('ECMWF ') 00524 CALL READ_GRIB(HGRIB,KLUOUT,170,IRET,PT2) 00525 CASE ('ARPEGE','ALADIN','MOCAGE') 00526 ILTYPE=111 00527 CALL READ_GRIB(HGRIB,KLUOUT,11,IRET,PT2,KLTYPE=ILTYPE) 00528 IF (IRET /= 0) THEN 00529 ILTYPE=105 00530 CALL READ_GRIB(HGRIB,KLUOUT,11,IRET,PT2,KLTYPE=ILTYPE) 00531 ENDIF 00532 CASE DEFAULT 00533 CALL ABOR1_SFX('MODE_READ_GRIB:READ_GRIB_T2:OPTION NOT SUPPORTED '//HINMODEL) 00534 END SELECT 00535 ! 00536 IF (IRET /= 0) THEN 00537 CALL ABOR1_SFX('MODE_READ_GRIB: DEEP SOIL TEMPERATURE MISSING (READ_GRIB_T2)') 00538 END IF 00539 ! 00540 IF (SIZE(PMASK)==SIZE(PT2)) & 00541 WHERE (PMASK(:)/=1.) PT2 = XUNDEF 00542 ! 00543 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_T2',1,ZHOOK_HANDLE) 00544 END SUBROUTINE READ_GRIB_T2 00545 !------------------------------------------------------------------- 00546 !------------------------------------------------------------------- 00547 SUBROUTINE PUT_LAYER_DEPTH(KLUOUT,KLEV,HROUT,KLTYPE,KLEV1,KLEV2,KNLAYERDEEP,PV4,PV,PD) 00548 ! 00549 IMPLICIT NONE 00550 ! 00551 INTEGER, INTENT(IN) :: KLUOUT 00552 INTEGER, INTENT(IN) :: KLEV 00553 CHARACTER(LEN=*), INTENT(IN) :: HROUT 00554 INTEGER, INTENT(INOUT) :: KLTYPE 00555 INTEGER, INTENT(IN) :: KLEV1 00556 INTEGER, INTENT(IN) :: KLEV2 00557 INTEGER, INTENT(IN) :: KNLAYERDEEP 00558 REAL, INTENT(IN) :: PV4 00559 REAL, INTENT(IN) :: PV 00560 REAL, INTENT(OUT) :: PD 00561 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00562 ! 00563 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:PUT_LAYER_DEPTH',0,ZHOOK_HANDLE) 00564 ! 00565 IF (KLEV2 == -1) KLTYPE = 0 00566 IF (KLTYPE==112) THEN 00567 PD = (KLEV2 - KLEV1) / 100. 00568 ELSE 00569 IF (KNLAYERDEEP == 4) THEN 00570 PD = PV4 00571 ELSE 00572 PD = PV 00573 END IF 00574 WRITE (KLUOUT,'(A,i1,A,f5.2,A)') 'MODE_READ_GRIB:'//TRIM(HROUT)//': | Level ',& 00575 KLEV,' height not available, assuming ',PD,' m' 00576 END IF 00577 ! 00578 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:PUT_LAYER_DEPTH',1,ZHOOK_HANDLE) 00579 END SUBROUTINE PUT_LAYER_DEPTH 00580 !------------------------------------------------------------------- 00581 ! ####################### 00582 SUBROUTINE FILL_PFIELD(KLUOUT,HROUT,KNLAYERDEEP,PDIN,PFIELDIN,PMASK,PFIELDOUT,PDOUT) 00583 ! ####################### 00584 ! 00585 USE MODD_SURF_PAR, ONLY : XUNDEF 00586 ! 00587 IMPLICIT NONE 00588 ! 00589 INTEGER, INTENT(IN) :: KLUOUT 00590 CHARACTER(LEN=*), INTENT(IN) :: HROUT 00591 INTEGER, INTENT(IN) :: KNLAYERDEEP 00592 REAL, DIMENSION(:), INTENT(IN) :: PDIN 00593 REAL, DIMENSION(:,:), INTENT(IN) :: PFIELDIN 00594 REAL, DIMENSION(:), INTENT(IN) :: PMASK 00595 REAL, DIMENSION(:,:), POINTER :: PFIELDOUT 00596 REAL, DIMENSION(:,:), POINTER :: PDOUT 00597 ! 00598 CHARACTER(LEN=20) :: FMT0 00599 INTEGER :: JL 00600 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00601 ! 00602 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:FILL_PFIELD',0,ZHOOK_HANDLE) 00603 !-------------------------------------------------------------------------------- 00604 ! 1. Display the number of layer found 00605 ! ----------------------- 00606 WRITE(FMT0,FMT='(A8,I1,A11)') '(A,I1,A,',KNLAYERDEEP,'(F5.2,","))' 00607 WRITE (KLUOUT,FMT=FMT0) 'MODE_READ_GRIB:'//TRIM(HROUT)//': | ',KNLAYERDEEP,& 00608 ' deep layers, heights are : ',PDIN(1:KNLAYERDEEP) 00609 !-------------------------------------------------------------------------------- 00610 ! 2. Set temperature profile and layer thicknesses 00611 ! ----------------------------------------------- 00612 ALLOCATE(PFIELDOUT(SIZE(PFIELDIN,1),KNLAYERDEEP)) 00613 ALLOCATE(PDOUT(SIZE(PFIELDIN,1),KNLAYERDEEP)) 00614 ! 00615 DO JL=1,KNLAYERDEEP 00616 PDOUT(:,JL)=SUM(PDIN(1:JL)) 00617 PFIELDOUT(:,JL)=PFIELDIN(:,JL) 00618 IF (SIZE(PMASK)==SIZE(PFIELDOUT,1)) & 00619 WHERE (PMASK(:)/=1.) PFIELDOUT(:,JL) = XUNDEF 00620 ENDDO 00621 ! 00622 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:FILL_PFIELD',1,ZHOOK_HANDLE) 00623 END SUBROUTINE FILL_PFIELD 00624 !------------------------------------------------------------------- 00625 ! ####################### 00626 SUBROUTINE READ_GRIB_TG_ECMWF(HGRIB,KLUOUT,HINMODEL,PMASK,PTG,PD) 00627 ! ####################### 00628 ! 00629 IMPLICIT NONE 00630 ! 00631 !* dummy arguments 00632 ! --------------- 00633 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00634 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00635 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00636 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 00637 REAL, DIMENSION(:,:), POINTER :: PTG ! field to initialize 00638 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 00639 ! 00640 !* local variables 00641 ! --------------- 00642 INTEGER(KIND=kindOfInt) :: IRET ! return code 00643 INTEGER :: ILTYPE ! type of level (Grib code table 3) 00644 INTEGER :: ILEV1 ! level definition 00645 INTEGER :: ILEV2 ! level definition 00646 INTEGER :: JL ! layer loop counter 00647 INTEGER :: INLAYERDEEP! number of deep moisture layers 00648 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() ! first layer temperature 00649 REAL, DIMENSION(:,:), ALLOCATABLE:: ZTG ! first layer temperature 00650 REAL, DIMENSION(:) , ALLOCATABLE:: ZD 00651 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00652 !-------------------------------------------------------------------------------- 00653 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TG_ECMWF',0,ZHOOK_HANDLE) 00654 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_TG_ECMWF: | Reading soil temperature' 00655 ! 00656 ALLOCATE(ZD(5)) 00657 ! 00658 ! 1. Search and read level 1 (and its depth) 00659 ! -------------------------------------- 00660 ILTYPE= -1 00661 ILEV1 = -1 00662 ILEV2 = -1 00663 CALL READ_GRIB(HGRIB,KLUOUT,139,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 00664 ! 00665 IF (IRET== 0) THEN 00666 CALL PUT_LAYER_DEPTH(KLUOUT,1,'READ_GRIB_TG_ECMWF',ILTYPE,ILEV1,ILEV2,4,0.07,0.07,ZD(1)) 00667 ALLOCATE(ZTG(SIZE(ZFIELD),5)) 00668 ZTG(:,1)=ZFIELD 00669 ELSE 00670 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL TEMPERATURE LEVEL 1 MISSING (READ_GRIB_TG_ECMWF)') 00671 ENDIF 00672 ! 00673 ! 2. Search and read level 4 (and its depth) This level is optionnal 00674 ! --------------------------------------------------------------- 00675 ILTYPE= -1 00676 ILEV1 = -1 00677 ILEV2 = -1 00678 CALL READ_GRIB(HGRIB,KLUOUT,236,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 00679 ! 00680 IF (IRET == 0) THEN 00681 INLAYERDEEP = 4 00682 CALL PUT_LAYER_DEPTH(KLUOUT,4,'READ_GRIB_TG_ECMWF',ILTYPE,ILEV1,ILEV2,INLAYERDEEP,1.89,1.89,ZD(4)) 00683 ZTG(:,4)=ZFIELD 00684 ELSE 00685 INLAYERDEEP = 3 00686 ZD(4) = 0. 00687 ENDIF 00688 ! 00689 ! 3. Search and read level 3 (and its depth) This level is optionnal 00690 ! --------------------------------------------------------------- 00691 ILTYPE= -1 00692 ILEV1 = -1 00693 ILEV2 = -1 00694 CALL READ_GRIB(HGRIB,KLUOUT,183,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 00695 ! 00696 IF (IRET == 0) THEN 00697 CALL PUT_LAYER_DEPTH(KLUOUT,3,'READ_GRIB_TG_ECMWF',ILTYPE,ILEV1,ILEV2,INLAYERDEEP,0.72,0.42,ZD(3)) 00698 ZTG(:,3)=ZFIELD 00699 ELSE 00700 INLAYERDEEP = 2 00701 ZD(3) = 0. 00702 ENDIF 00703 ! 00704 ! 4. Search and read level 2 (and its depth) 00705 ! --------------------------------------- 00706 ILTYPE= -1 00707 ILEV1 = -1 00708 ILEV2 = -1 00709 CALL READ_GRIB(HGRIB,KLUOUT,170,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 00710 ! 00711 IF (IRET== 0) THEN 00712 CALL PUT_LAYER_DEPTH(KLUOUT,2,'READ_GRIB_TG_ECMWF',ILTYPE,ILEV1,ILEV2,INLAYERDEEP,0.21,0.42,ZD(2)) 00713 ZTG(:,2)=ZFIELD 00714 DEALLOCATE(ZFIELD) 00715 ELSE 00716 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL TEMPERATURE LEVEL 2 MISSING (READ_GRIB_TG_ECMWF)') 00717 ENDIF 00718 !-------------------------------------------------------------------------------- 00719 ! 5. Assumes uniform temperature profile up to 3m depth 00720 ! ------------------------------------------------- 00721 ! 00722 IF(SUM(ZD(1:INLAYERDEEP)) < 3.) THEN 00723 !We add a temperature layer 00724 INLAYERDEEP=INLAYERDEEP+1 00725 ZD(INLAYERDEEP)=3.-SUM(ZD(1:INLAYERDEEP-1)) 00726 ZTG(:,INLAYERDEEP)=ZTG(:,INLAYERDEEP-1) 00727 ENDIF 00728 ! 00729 !-------------------------------------------------------------------------------- 00730 ! 6. Set temperature profile and layer thicknesses 00731 ! ---------------------------------------------- 00732 CALL FILL_PFIELD(KLUOUT,'READ_GRIB_TG_ECMWF',INLAYERDEEP,ZD,ZTG,PMASK,PTG,PD) 00733 DEALLOCATE(ZD) 00734 DEALLOCATE(ZTG) 00735 ! 00736 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TG_ECMWF',1,ZHOOK_HANDLE) 00737 END SUBROUTINE READ_GRIB_TG_ECMWF 00738 !------------------------------------------------------------------- 00739 ! ####################### 00740 SUBROUTINE READ_GRIB_WG_ECMWF_1(HGRIB,KLUOUT,HINMODEL,PMASK,PWG,PD) 00741 ! ####################### 00742 ! 00743 ! This tasks is divided in the following steps : 00744 ! - computing the MesoNH constants 00745 ! - reading the grib datas according to the type of file (ECMWF/Arpege/Aladin) 00746 ! - converting from specific humidity to relative humidity 00747 ! - interpolation with land mask 00748 ! - converting back from relative humidity to specific humidity with MesoNH constants 00749 ! Five different models are supported : 00750 ! - ECMWF with 2 layers (untested) 00751 ! - ECMWF with 3 layers (archive before 1991 - Blondin model) 00752 ! - ECMWF with 4 layers (archive after 1991 - Viterbo model) 00753 ! - Arpege/Aladin before ISBA (I don't know the name of this model) 00754 ! - Arpege/Aladin with ISBA model 00755 ! The available model is detect according to the fields which are presents : 00756 ! - ECMWF archive : loads as many layers as possible 00757 ! - Arpege/Aladin archive : ISBA model needs Clay and Sans fraction fields, if they 00758 ! are present, they are used and the model is declared to be ISBA. 00759 ! To detect the height of the layers, two methods are used : 00760 ! - if level type is not 112, a default value is assumed and a warning message is 00761 ! displayed 00762 ! - if level type is ID 112, then the position of the top and bottom surface may be 00763 ! given. If they are present, they are used, if not the default value is taken and 00764 ! a warning message is issued. 00765 ! 00766 IMPLICIT NONE 00767 ! 00768 !* dummy arguments 00769 ! --------------- 00770 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00771 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00772 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00773 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 00774 REAL, DIMENSION(:,:), POINTER :: PWG ! field to initialize 00775 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 00776 ! 00777 !* local variables 00778 ! --------------- 00779 INTEGER(KIND=kindOfInt) :: IRET ! return code 00780 INTEGER :: IPAR ! parameter number for field reading 00781 INTEGER :: ILTYPE ! type of level (Grib code table 3) 00782 INTEGER :: ILEV1 ! level definition 00783 INTEGER :: ILEV2 ! level definition 00784 INTEGER :: INLAYERDEEP! number of deep moisture layers 00785 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() ! first layer temperature 00786 REAL, DIMENSION(:,:), ALLOCATABLE:: ZWG ! first layer temperature 00787 REAL, DIMENSION(:) , ALLOCATABLE:: ZD 00788 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00789 !-------------------------------------------------------------------------------- 00790 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WG_ECMWF_1',0,ZHOOK_HANDLE) 00791 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_WG_ECMWF_1: | Reading soil moisture' 00792 ! 00793 ALLOCATE(ZD(4)) 00794 ! 00795 ! 1. Search and read level 1 (and its depth) 00796 ! -------------------------------------- 00797 ILTYPE= -1 00798 ILEV1 = -1 00799 ILEV2 = -1 00800 IPAR=39 00801 CALL READ_GRIB(HGRIB,KLUOUT,140,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2,KPARAM2=IPAR) 00802 ! 00803 IF (IRET == 0) THEN 00804 CALL PUT_LAYER_DEPTH(KLUOUT,1,'READ_GRIB_WG_ECMWF_1',ILTYPE,ILEV1,ILEV2,4,0.07,0.07,ZD(1)) 00805 ALLOCATE(ZWG(SIZE(ZFIELD,1),4)) 00806 ZWG(:,1)=ZFIELD 00807 ! 00808 IF (IPAR==140) ZWG(:,1)=ZWG(:,1) / ZD(1) 00809 ELSE 00810 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL MOISTURE LEVEL 1 MISSING (READ_GRIB_WG_ECMWF_1)') 00811 ENDIF 00812 ! 00813 ! 2. Search and read level 4 (and its depth) This level is optionnal 00814 ! --------------------------------------------------------------- 00815 ILTYPE= -1 00816 ILEV1 = -1 00817 ILEV2 = -1 00818 IPAR=42 00819 CALL READ_GRIB(HGRIB,KLUOUT,237,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2,KPARAM2=IPAR) 00820 ! 00821 IF (IRET == 0) THEN 00822 INLAYERDEEP = 4 00823 CALL PUT_LAYER_DEPTH(KLUOUT,4,'READ_GRIB_WG_ECMWF_1',ILTYPE,ILEV1,ILEV2,INLAYERDEEP,1.89,1.89,ZD(4)) 00824 ZWG(:,4)=ZFIELD 00825 ! 00826 IF (IPAR==237) ZWG(:,4)=ZWG(:,4) / ZD(1) 00827 ELSE 00828 INLAYERDEEP = 3 00829 ZD(4) = 0. 00830 ENDIF 00831 ! 00832 ! 3. Search and read level 3 (and its depth) This level is optionnal 00833 ! --------------------------------------------------------------- 00834 ILTYPE= -1 00835 ILEV1 = -1 00836 ILEV2 = -1 00837 IPAR=41 00838 CALL READ_GRIB(HGRIB,KLUOUT,184,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2,KPARAM2=IPAR) 00839 ! 00840 IF (IRET == 0) THEN 00841 CALL PUT_LAYER_DEPTH(KLUOUT,3,'READ_GRIB_WG_ECMWF_1',ILTYPE,ILEV1,ILEV2,INLAYERDEEP,0.72,0.42,ZD(3)) 00842 ZWG(:,3)=ZFIELD 00843 ! 00844 IF (IPAR==184) ZWG(:,3)=ZWG(:,3) / ZD(1) 00845 ELSE 00846 INLAYERDEEP = 2 00847 ZD(3) = 0. 00848 ENDIF 00849 ! 00850 ! 4. Search and read level 2 (and its depth) 00851 ! --------------------------------------- 00852 ILTYPE= -1 00853 ILEV1 = -1 00854 ILEV2 = -1 00855 IPAR=40 00856 CALL READ_GRIB(HGRIB,KLUOUT,171,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2,KPARAM2=IPAR) 00857 ! 00858 IF (IRET == 0) THEN 00859 CALL PUT_LAYER_DEPTH(KLUOUT,2,'READ_GRIB_WG_ECMWF_1',ILTYPE,ILEV1,ILEV2,INLAYERDEEP,0.21,0.42,ZD(2)) 00860 ZWG(:,2)=ZFIELD 00861 DEALLOCATE(ZFIELD) 00862 ! 00863 IF (IPAR==171) ZWG(:,2)=ZWG(:,2) / ZD(1) 00864 ELSE 00865 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL MOISTURE LEVEL 2 MISSING (READ_GRIB_WG_ECMWF_1)') 00866 ENDIF 00867 ! 00868 !-------------------------------------------------------------------------------- 00869 ! 00870 CALL FILL_PFIELD(KLUOUT,'READ_GRIB_WG_ECMWF_1',INLAYERDEEP,ZD,ZWG,PMASK,PWG,PD) 00871 DEALLOCATE(ZD) 00872 DEALLOCATE(ZWG) 00873 ! 00874 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WG_ECMWF_1',1,ZHOOK_HANDLE) 00875 END SUBROUTINE READ_GRIB_WG_ECMWF_1 00876 !------------------------------------------------------------------- 00877 ! ####################### 00878 SUBROUTINE ECMWF_WGI(PTG,PWG,PWGI) 00879 ! ####################### 00880 ! 00881 ! ECMWF grib only contain (ice+water) content. 00882 ! This routine computes iced part and water part according to the formula 00883 ! given in ECMWF documentation. But we use real water content instead of 00884 ! (CL+CH) times saturation capacity. 00885 ! 00886 USE MODD_CSTS, ONLY : XTT, XPI 00887 ! 00888 IMPLICIT NONE 00889 ! 00890 !* dummy arguments 00891 ! --------------- 00892 REAL, DIMENSION(:,:), INTENT(IN) :: PTG ! Temperature profil 00893 REAL, DIMENSION(:,:), INTENT(INOUT) :: PWG ! INPUT contains (water+ice) profil 00894 ! OUTPUT contains only water profil 00895 REAL, DIMENSION(:,:), INTENT(OUT) :: PWGI ! ice profil 00896 ! 00897 !* local variables 00898 ! --------------- 00899 REAL :: ZT1, ZT2 ! Temperature threshold 00900 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00901 !-------------------------------------------------------------------------------- 00902 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:ECMWF_WGI',0,ZHOOK_HANDLE) 00903 ! 00904 ZT1=XTT + 1. 00905 ZT2=XTT - 3. 00906 ! 00907 WHERE(PTG(:,:) > ZT1) 00908 PWGI(:,:) = 0. 00909 ELSEWHERE(PTG(:,:) < ZT2) 00910 PWGI(:,:) = PWG(:,:) 00911 PWG(:,:) = 0. 00912 ELSEWHERE 00913 PWGI(:,:)=PWG(:,:) * 0.5* (1 - sin(XPI * (PTG(:,:) - 0.5*ZT1 - 0.5*ZT2) / & 00914 (ZT1 - ZT2 ) )) 00915 PWG(:,:) = PWG(:,:) - PWGI(:,:) 00916 ENDWHERE 00917 ! 00918 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:ECMWF_WGI',1,ZHOOK_HANDLE) 00919 END SUBROUTINE ECMWF_WGI 00920 !-------------------------------------------------------------------------------- 00921 ! ####################### 00922 SUBROUTINE HARMONIZE_GRIB_WG_WGI_ECMWF(HGRIB,KLUOUT,HINMODEL,PMASK,PWG,PD,PWGI) 00923 ! ####################### 00924 ! 00925 ! ECMWF grib only contain (ice+water) content. 00926 ! This routine computes iced part and water part according to the formula 00927 ! given in ECMWF documentation. But we use real water content instead of 00928 ! (CL+CH) times saturation capacity. 00929 ! 00930 USE MODD_CSTS, ONLY : XTT, XPI 00931 ! 00932 IMPLICIT NONE 00933 ! 00934 !* dummy arguments 00935 ! --------------- 00936 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 00937 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 00938 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 00939 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 00940 REAL, DIMENSION(:,:), OPTIONAL, POINTER :: PWG ! INPUT contains (water+ice) profil 00941 ! OUTPUT contains only water profil 00942 REAL, DIMENSION(:,:), OPTIONAL, POINTER :: PD ! thickness of each layer 00943 REAL, DIMENSION(:,:), OPTIONAL, POINTER :: PWGI ! ice profil 00944 !* local variables 00945 ! --------------- 00946 REAL, DIMENSION(:,:), POINTER :: ZWG => NULL() ! profile of soil water contents 00947 REAL, DIMENSION(:,:), POINTER :: ZD => NULL() ! thickness of each layer 00948 REAL, DIMENSION(:,:), POINTER :: ZTG => NULL() ! profile of temperature 00949 REAL, DIMENSION(:,:), POINTER :: ZDT => NULL() ! thickness of each temperature layer 00950 REAL, DIMENSION(:,:), ALLOCATABLE:: ZWGI ! profile of soil ice contents 00951 REAL :: ZT1, ZT2 ! Temperature threshold 00952 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00953 !-------------------------------------------------------------------------------- 00954 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:HARMONIZE_GRIB_WG_WGI_ECMWF',0,ZHOOK_HANDLE) 00955 ! 00956 CALL READ_GRIB_TG_ECMWF(HGRIB,KLUOUT,HINMODEL,PMASK,ZTG,ZDT) 00957 CALL READ_GRIB_WG_ECMWF_1(HGRIB,KLUOUT,HINMODEL,PMASK,ZWG,ZD) 00958 ! 00959 IF (SIZE(ZTG,2) .LT. SIZE(ZWG,2)) THEN 00960 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:HARMONIZE_GRIB_WG_WGI_ECMWF: ' 00961 WRITE (KLUOUT,'(A)') 'ERROR, YOU HAVE NOT THE SAME NUMBER OF LEVELS ' 00962 WRITE (KLUOUT,'(A)') 'IN SOIL FOR TEMPERATURE AND HUMIDITY ' 00963 WRITE (KLUOUT,'(A)') 'VERIFY GRIB FILE ' 00964 CALL ABOR1_SFX("MODE_READ_GRIB:HARMONIZE_GRIB_WG_WGI_ECMWF: VERIFY NUMBER OF LEVELS IN GRIB FILE") 00965 ENDIF 00966 ! 00967 IF (PRESENT(PD)) THEN 00968 ALLOCATE(PD(SIZE(ZD,1),SIZE(ZD,2))) 00969 PD(:,:)=ZD(:,:) 00970 ENDIF 00971 IF (PRESENT(PWGI)) THEN 00972 ALLOCATE(PWGI(SIZE(ZWG,1),SIZE(ZWG,2))) 00973 PWGI(:,:)=0. 00974 ENDIF 00975 ! 00976 !If same vertical grids are taken into account for WG and TG we can 00977 !compute ice content and new water content 00978 IF(ALL(ZDT(:,1:SIZE(ZWG,2))==ZD(:,1:SIZE(ZWG,2)))) THEN 00979 ALLOCATE(ZWGI(SIZE(ZWG,1),SIZE(ZWG,2))) 00980 CALL ECMWF_WGI(ZTG(:,1:SIZE(ZWG,2)),ZWG,ZWGI) 00981 IF (PRESENT(PWGI)) PWGI(:,:)=ZWGI(:,:) 00982 DEALLOCATE(ZWGI) 00983 ENDIF 00984 ! 00985 IF (PRESENT(PWG)) THEN 00986 ALLOCATE(PWG(SIZE(ZWG,1),SIZE(ZWG,2))) 00987 PWG(:,:)=ZWG(:,:) 00988 ENDIF 00989 ! 00990 DEALLOCATE(ZWG) 00991 DEALLOCATE(ZD) 00992 DEALLOCATE(ZTG) 00993 DEALLOCATE(ZDT) 00994 ! 00995 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:HARMONIZE_GRIB_WG_WGI_ECMWF',1,ZHOOK_HANDLE) 00996 END SUBROUTINE HARMONIZE_GRIB_WG_WGI_ECMWF 00997 !-------------------------------------------------------------------------------- 00998 ! ####################### 00999 SUBROUTINE READ_GRIB_WG_ECMWF(HGRIB,KLUOUT,HINMODEL,PMASK,PFIELD,PD) 01000 ! ####################### 01001 ! 01002 USE MODD_GRID_GRIB, ONLY : NNI 01003 USE MODD_SURF_PAR, ONLY : XUNDEF 01004 ! 01005 IMPLICIT NONE 01006 ! 01007 !* dummy arguments 01008 ! --------------- 01009 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01010 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01011 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01012 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01013 REAL, DIMENSION(:,:), POINTER :: PFIELD ! field to initialize 01014 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 01015 !* local variables 01016 ! --------------- 01017 INTEGER(KIND=kindOfInt) :: IRET ! return code 01018 REAL, DIMENSION(:), POINTER :: ZSLT => NULL() ! soil type 01019 REAL, DIMENSION(:), ALLOCATABLE:: ZWWILT ! ECMWF wilting point 01020 REAL, DIMENSION(:), ALLOCATABLE:: ZWFC ! ECMWF field capacity 01021 INTEGER :: JL ! loop counter on layers 01022 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01023 ! 01024 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WG_ECMWF',0,ZHOOK_HANDLE) 01025 ! 01026 CALL HARMONIZE_GRIB_WG_WGI_ECMWF(HGRIB,KLUOUT,HINMODEL,PMASK,PWG=PFIELD,PD=PD) 01027 ! 01028 ! 1. Get soil type to compute SWI 01029 ! ---------------------------- 01030 CALL READ_GRIB(HGRIB,KLUOUT,43,IRET,ZSLT) 01031 !-------------------------------------------------------------------------------- 01032 ALLOCATE (ZWFC(SIZE(PFIELD,1))) 01033 ALLOCATE (ZWWILT(SIZE(PFIELD,1))) 01034 ZWFC (:) = 0. 01035 ZWWILT(:) = 0. 01036 ! 01037 IF (IRET == 0) THEN 01038 ! 01039 ! 2.1 Convert from specific humidity to relative humidity using soil types 01040 ! -------------------------------------------------------------------- 01041 WHERE (ZSLT(:)==1.) 01042 ZWFC(:) = 0.242 01043 ZWWILT(:) = 0.059 01044 ELSEWHERE (ZSLT(:)==2.) 01045 ZWFC(:) = 0.346 01046 ZWWILT(:) = 0.151 01047 ELSEWHERE (ZSLT(:)==3.) 01048 ZWFC(:) = 0.382 01049 ZWWILT(:) = 0.133 01050 ELSEWHERE (ZSLT(:)==4.) 01051 ZWFC(:) = 0.448 01052 ZWWILT(:) = 0.279 01053 ELSEWHERE (ZSLT(:)==5.) 01054 ZWFC(:) = 0.541 01055 ZWWILT(:) = 0.335 01056 ELSEWHERE (ZSLT(:)==6.) 01057 ZWFC(:) = 0.662 01058 ZWWILT(:) = 0.267 01059 ENDWHERE 01060 ! 01061 ELSE 01062 ! 01063 ! 2.2 Convert from specific humidity to relative humidity single soil type 01064 ! -------------------------------------------------------------------- 01065 ! Compute model's constants 01066 IF (SIZE(PFIELD,2)==4) THEN 01067 ZWFC(:) = 0.323 01068 ZWWILT(:) = 0.171 01069 ELSE 01070 ZWFC(:) = 0.171 01071 ZWWILT(:) = 0.086 01072 END IF 01073 ! 01074 ENDIF 01075 ! 01076 DO JL=1,SIZE(PFIELD,2) 01077 WHERE ( PFIELD(:,JL).NE.XUNDEF .AND. ZWFC(:).NE.0. ) 01078 PFIELD(:,JL) = (PFIELD(:,JL) - ZWWILT(:)) / (ZWFC(:) - ZWWILT(:)) 01079 ELSEWHERE 01080 PFIELD(:,JL) = 0. 01081 ENDWHERE 01082 ENDDO 01083 ! 01084 IF (ASSOCIATED(ZSLT)) DEALLOCATE(ZSLT) 01085 DEALLOCATE(ZWFC) 01086 DEALLOCATE(ZWWILT) 01087 ! 01088 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WG_ECMWF',1,ZHOOK_HANDLE) 01089 END SUBROUTINE READ_GRIB_WG_ECMWF 01090 !---------------------------------------------------------------------------- 01091 ! ####################### 01092 SUBROUTINE READ_GRIB_WGI_ECMWF(HGRIB,KLUOUT,HINMODEL,PMASK,PFIELD,PD) 01093 ! ####################### 01094 ! 01095 USE MODD_GRID_GRIB, ONLY : NNI 01096 USE MODD_SURF_PAR, ONLY : XUNDEF 01097 ! 01098 IMPLICIT NONE 01099 ! 01100 !* dummy arguments 01101 ! --------------- 01102 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01103 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01104 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01105 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01106 REAL, DIMENSION(:,:), POINTER :: PFIELD ! field to initialize 01107 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 01108 !* local variables 01109 ! --------------- 01110 INTEGER(KIND=kindOfInt) :: IRET ! return code 01111 REAL, DIMENSION(:), POINTER :: ZSLT => NULL() ! soil type 01112 REAL, DIMENSION(:) , ALLOCATABLE:: ZWSAT ! ECMWF saturation 01113 INTEGER :: JL ! loop counter on layers 01114 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01115 !-------------------------------------------------------------------------------- 01116 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WGI_ECMWF',0,ZHOOK_HANDLE) 01117 ! 01118 CALL HARMONIZE_GRIB_WG_WGI_ECMWF(HGRIB,KLUOUT,HINMODEL,PMASK,PD=PD,PWGI=PFIELD) 01119 ! 01120 ! 1. Get soil type to compute WSAT 01121 ! ---------------------------- 01122 CALL READ_GRIB(HGRIB,KLUOUT,43,IRET,ZSLT) 01123 !-------------------------------------------------------------------------------- 01124 ALLOCATE (ZWSAT(SIZE(PFIELD,1))) 01125 ZWSAT(:)=0. 01126 ! 01127 IF (IRET == 0) THEN 01128 ! 01129 ! 2.1 Convert from specific humidity to relative humidity using soil types 01130 ! -------------------------------------------------------------------- 01131 WHERE (ZSLT(:)==1.) 01132 ZWSAT(:) = 0.403 01133 ELSEWHERE (ZSLT(:)==2.) 01134 ZWSAT(:) = 0.439 01135 ELSEWHERE (ZSLT(:)==3.) 01136 ZWSAT(:) = 0.430 01137 ELSEWHERE (ZSLT(:)==4.) 01138 ZWSAT(:) = 0.520 01139 ELSEWHERE (ZSLT(:)==5.) 01140 ZWSAT(:) = 0.614 01141 ELSEWHERE (ZSLT(:)==6.) 01142 ZWSAT(:) = 0.766 01143 ENDWHERE 01144 ! 01145 ELSE 01146 ! 01147 ! 2.2 Convert from specific humidity to relative humidity single soil type 01148 ! -------------------------------------------------------------------- 01149 ! Compute model's constants 01150 IF (SIZE(PFIELD,2)==4) THEN 01151 ZWSAT(:) = 0.472 01152 ELSE 01153 ZWSAT(:) = 0.286 01154 END IF 01155 ! 01156 ENDIF 01157 ! 01158 ! Then perform conversion 01159 DO JL=1,SIZE(PFIELD,2) 01160 WHERE ( PFIELD(:,JL).NE.XUNDEF .AND. ZWSAT(:).NE.0. ) 01161 PFIELD(:,JL) = PFIELD(:,JL) / ZWSAT(:) 01162 ELSEWHERE 01163 PFIELD(:,JL) = 0. 01164 ENDWHERE 01165 ENDDO 01166 ! 01167 IF (ASSOCIATED(ZSLT)) DEALLOCATE(ZSLT) 01168 DEALLOCATE(ZWSAT) 01169 ! 01170 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WGI_ECMWF',1,ZHOOK_HANDLE) 01171 END SUBROUTINE READ_GRIB_WGI_ECMWF 01172 !------------------------------------------------------------------- 01173 !------------------------------------------------------------------- 01174 ! ####################### 01175 SUBROUTINE READ_GRIB_TG_METEO_FRANCE(HGRIB,KLUOUT,HINMODEL,PMASK,PTG,PDT) 01176 ! ####################### 01177 ! 01178 USE MODD_SURF_PAR, ONLY : XUNDEF 01179 ! 01180 IMPLICIT NONE 01181 ! 01182 !* dummy arguments 01183 ! --------------- 01184 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01185 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01186 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01187 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01188 REAL, DIMENSION(:,:), POINTER :: PTG ! field to initialize 01189 REAL, DIMENSION(:,:), POINTER :: PDT ! thickness of each layer 01190 !* local variables 01191 ! --------------- 01192 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() ! field to read 01193 INTEGER :: JL ! layer loop counter 01194 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01195 !-------------------------------------------------------------------------------- 01196 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TG_METEO_FRANCE',0,ZHOOK_HANDLE) 01197 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_TG_METEO_FRANCE: | Reading soil temperature' 01198 !-------------------------------------------------------------------------------- 01199 ! 1. Allocate soil temperature profile 01200 ! --------------------------------- 01201 !-------------------------------------------------------------------------------- 01202 ! 2. Search and read level 1 (and its depth) 01203 ! --------------------------------------- 01204 CALL READ_GRIB_TS(HGRIB,KLUOUT,HINMODEL,PMASK,ZFIELD) 01205 ! 01206 ALLOCATE(PTG(SIZE(ZFIELD),3)) 01207 ALLOCATE(PDT(SIZE(ZFIELD),3)) 01208 ! 01209 PTG(:,1) = ZFIELD(:) 01210 PDT(:,1) = 0. 01211 !-------------------------------------------------------------------------------- 01212 ! 3. Deep soil temperature 01213 ! --------------------- 01214 CALL READ_GRIB_T2(HGRIB,KLUOUT,HINMODEL,PMASK,ZFIELD) 01215 ! 01216 PTG(:,2) = ZFIELD(:) 01217 PDT(:,2) = 0.2 ! deep temperature depth assumed equal to 0.2m 01218 DEALLOCATE(ZFIELD) 01219 !-------------------------------------------------------------------------------- 01220 ! 4. Assumes uniform temperature profile below 01221 ! ----------------------------------------- 01222 PTG(:,3) = PTG(:,2) 01223 PDT(:,3) = 3. ! temperature profile down to 3m 01224 ! 01225 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TG_METEO_FRANCE',1,ZHOOK_HANDLE) 01226 END SUBROUTINE READ_GRIB_TG_METEO_FRANCE 01227 !------------------------------------------------------------------- 01228 ! ####################### 01229 SUBROUTINE READ_GRIB_SAND_CLAY_METEO_FRANCE(HGRIB,KLUOUT,HINMODEL,PSAND,PCLAY,GISBA) 01230 ! ###################### 01231 ! 01232 IMPLICIT NONE 01233 ! 01234 !* dummy arguments 01235 ! --------------- 01236 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01237 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01238 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01239 REAL, DIMENSION(:), POINTER :: PSAND ! field to initialize 01240 REAL, DIMENSION(:), POINTER :: PCLAY ! thickness of each layer 01241 LOGICAL, INTENT(OUT) :: GISBA ! T: surface scheme in file is ISBA 01242 !* local variables 01243 ! --------------- 01244 INTEGER(KIND=kindOfInt) :: IRET ! return code 01245 INTEGER :: IPAR ! parameter number for field reading 01246 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01247 !------------------------------------------------------------------------------- 01248 ! 1. Search and read clay fraction if available 01249 ! ------------------------------------------ 01250 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_SAND_CLAY_METEO_FRANCE',0,ZHOOK_HANDLE) 01251 ! 01252 IF (HINMODEL == 'ARPEGE' .OR. HINMODEL == 'MOCAGE') IPAR=171 01253 IF (HINMODEL == 'ALADIN') IPAR=128 01254 CALL READ_GRIB(HGRIB,KLUOUT,IPAR,IRET,PCLAY) 01255 ! 01256 ! if not available, the model is not ISBA (IWMODE=1) 01257 IF (IRET /= 0) THEN 01258 GISBA = .FALSE. 01259 ELSE 01260 GISBA = .TRUE. 01261 PCLAY(:) = PCLAY(:) / 100. ! this field is given in percent 01262 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_SAND_CLAY_METEO_FRANCE: | The soil model is ISBA' 01263 END IF 01264 !------------------------------------------------------------------------------- 01265 ! 2. Search and read sand fraction if available 01266 ! ------------------------------------------ 01267 IF (HINMODEL == 'ARPEGE' .OR. HINMODEL == 'MOCAGE') IPAR=172 01268 IF (HINMODEL == 'ALADIN') IPAR=129 01269 CALL READ_GRIB(HGRIB,KLUOUT,IPAR,IRET,PSAND) 01270 ! 01271 ! if not available, the model is not ISBA (IWMODE=1) 01272 IF (GISBA) THEN 01273 IF (IRET /= 0) THEN 01274 CALL ABOR1_SFX('MODE_READ_GRIB: SAND FRACTION MISSING (READ_GRIB_SAND_CLAY_METEO_FRANCE)') 01275 ELSE 01276 PSAND(:) = PSAND(:) / 100. ! this field is given in percent 01277 END IF 01278 END IF 01279 ! 01280 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_SAND_CLAY_METEO_FRANCE',1,ZHOOK_HANDLE) 01281 END SUBROUTINE READ_GRIB_SAND_CLAY_METEO_FRANCE 01282 !----------------------------------------------------------------------------------- 01283 ! ####################### 01284 SUBROUTINE READ_GRIB_WG_METEO_FRANCE(HGRIB,KLUOUT,HINMODEL,PMASK,PFIELD,PD) 01285 ! ####################### 01286 ! 01287 USE MODD_GRID_GRIB, ONLY : NNI 01288 USE MODD_SURF_PAR, ONLY : XUNDEF 01289 ! 01290 IMPLICIT NONE 01291 ! 01292 !* dummy arguments 01293 ! --------------- 01294 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01295 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01296 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01297 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01298 REAL, DIMENSION(:,:), POINTER :: PFIELD ! field to initialize 01299 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 01300 ! 01301 !* local variables 01302 ! --------------- 01303 LOGICAL :: GISBA ! T: surface scheme in file is ISBA 01304 INTEGER(KIND=kindOfInt) :: IRET ! return code 01305 INTEGER :: ILTYPE ! type of level (Grib code table 3) 01306 INTEGER :: ILEV1 ! level definition 01307 INTEGER :: ILEV2 ! level definition 01308 REAL, DIMENSION(:), POINTER :: ZCLAY => NULL() ! clay fraction 01309 REAL, DIMENSION(:), POINTER :: ZSAND => NULL() ! sand fraction 01310 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() 01311 REAL, DIMENSION(:), ALLOCATABLE:: ZWWILT ! wilting point 01312 REAL, DIMENSION(:), ALLOCATABLE:: ZWFC ! field capacity 01313 REAL, DIMENSION(:), ALLOCATABLE:: ZWSAT ! saturation 01314 INTEGER :: JL ! layer loop counter 01315 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01316 !------------------------------------------------------------------------------- 01317 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WG_METEO_FRANCE',0,ZHOOK_HANDLE) 01318 ! 01319 ! 1. Search and read clay and sand fractions if available 01320 ! ---------------------------------------------------- 01321 CALL READ_GRIB_SAND_CLAY_METEO_FRANCE(HGRIB,KLUOUT,HINMODEL,ZSAND,ZCLAY,GISBA) 01322 !------------------------------------------------------------------------------- 01323 IF (GISBA) THEN 01324 ALLOCATE(PFIELD(SIZE(ZSAND),3)) 01325 ALLOCATE(PD(SIZE(ZSAND),3)) 01326 ELSE 01327 ALLOCATE(PFIELD(NNI,3)) 01328 ALLOCATE(PD(NNI,3)) 01329 ENDIF 01330 ! 01331 PD(:,1) = 0. 01332 PD(:,2) = 0.20 01333 !------------------------------------------------------------------------------- 01334 ! 2. Read layer 1 moisture 01335 ! --------------------- 01336 ILEV1 = 0 01337 IF (HINMODEL == 'ARPEGE' .OR. HINMODEL=='MOCAGE') THEN 01338 ILTYPE = 112 01339 ILEV2 = 1 01340 CALL READ_GRIB(HGRIB,KLUOUT,153,IRET,ZFIELD,KLEV1=ILEV1,KLEV2=ILEV2) 01341 ELSE 01342 ILTYPE = 105 01343 ILEV2 = 0 01344 CALL READ_GRIB(HGRIB,KLUOUT,86,IRET,ZFIELD,KLEV1=ILEV1,KLEV2=ILEV2) 01345 ENDIF 01346 IF (IRET /= 0) THEN 01347 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL MOISTURE LEVEL 1 MISSING (READ_GRIB_WG_METEO_FRANCE)') 01348 END IF 01349 ! 01350 PFIELD(:,1) = ZFIELD(:) 01351 !------------------------------------------------------------------------------- 01352 ! 3. Read layer 2 moisture 01353 ! --------------------- 01354 IF (HINMODEL == 'ARPEGE' .OR. HINMODEL=='MOCAGE') THEN 01355 ILTYPE = 112 01356 ILEV1 = 0 01357 ILEV2 = 250 01358 CALL READ_GRIB(HGRIB,KLUOUT,153,IRET,ZFIELD,KLEV1=ILEV1,KLEV2=ILEV2) 01359 ELSE 01360 ILTYPE = 111 01361 ILEV1 = -1 01362 ILEV2 = -1 01363 CALL READ_GRIB(HGRIB,KLUOUT,86,IRET,ZFIELD,KLEV1=ILEV1,KLEV2=ILEV2) 01364 ENDIF 01365 IF (IRET /= 0) THEN 01366 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL MOISTURE LEVEL 2 MISSING (READ_GRIB_WG_METEO_FRANCE)') 01367 END IF 01368 ! 01369 PFIELD(:,2) = ZFIELD(:) 01370 !------------------------------------------------------------------------------- 01371 ! 4. Read layer 2 depth (ISBA only) 01372 ! ----------------------------- 01373 !* note that soil water reservoir is considered uniform between 0.2m and GRIB soil depth 01374 IF (GISBA) THEN 01375 IF (HINMODEL == 'ARPEGE' .OR. HINMODEL == 'MOCAGE') THEN 01376 CALL READ_GRIB(HGRIB,KLUOUT,173,IRET,ZFIELD) 01377 ELSE 01378 CALL READ_GRIB(HGRIB,KLUOUT,130,IRET,ZFIELD) 01379 ENDIF 01380 IF (IRET /= 0) THEN 01381 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL MOISTURE LEVEL 2 DEPTH MISSING (READ_GRIB_WG_METEO_FRANCE)') 01382 END IF 01383 PD(:,3) = ZFIELD(:) 01384 DEALLOCATE(ZFIELD) 01385 ELSE 01386 PD(:,3) = 2. 01387 END IF 01388 !------------------------------------------------------------------------------- 01389 ! 5. Compute relative humidity from units kg/m^2 01390 ! ------------------------------------------- 01391 ! Compute ISBA model constants (if needed) 01392 IF (GISBA) THEN 01393 ! 01394 !* updates Wg in m3/m3 01395 PFIELD(:,1) = PFIELD(:,1) / 10. 01396 PFIELD(:,2) = PFIELD(:,2) /(1000. * PD(:,3)) 01397 ! 01398 ALLOCATE (ZWSAT (SIZE(ZSAND))) 01399 ZWSAT (:) = (-1.08*100. * ZSAND(:) + 494.305) * 0.001 01400 PFIELD(:,1) = MAX(MIN(PFIELD(:,1),ZWSAT(:)),0.) 01401 PFIELD(:,2) = MAX(MIN(PFIELD(:,2),ZWSAT(:)),0.) 01402 DEALLOCATE(ZWSAT) 01403 DEALLOCATE (ZSAND) 01404 ! 01405 ALLOCATE (ZWWILT(SIZE(ZCLAY))) 01406 ALLOCATE (ZWFC (SIZE(ZCLAY))) 01407 ZWWILT(:) = 37.1342E-3 * SQRT( 100. * ZCLAY(:) ) 01408 ZWFC (:) = 89.0467E-3 * (100. * ZCLAY(:) )**0.3496 01409 PFIELD(:,1) = (PFIELD(:,1) - ZWWILT(:)) / (ZWFC(:) - ZWWILT(:)) 01410 PFIELD(:,2) = (PFIELD(:,2) - ZWWILT(:)) / (ZWFC(:) - ZWWILT(:)) 01411 DEALLOCATE (ZWWILT) 01412 DEALLOCATE (ZWFC) 01413 DEALLOCATE (ZCLAY) 01414 ! 01415 ELSE ! Non ISBA 01416 ! 01417 PFIELD(:,2) = (PFIELD(:,1)+PFIELD(:,2)) / (20. + 100.) 01418 PFIELD(:,1) = PFIELD(:,1) / 20. 01419 ! 01420 END IF 01421 ! 01422 PFIELD(:,3) = PFIELD(:,2) 01423 !-------------------------------------------------------------------------------- 01424 ! 6. Apply land mask 01425 ! --------------- 01426 IF (SIZE(PMASK)==SIZE(PFIELD,1)) THEN 01427 DO JL=1,SIZE(PFIELD,2) 01428 WHERE (PMASK(:)/=1.) PFIELD(:,JL) = XUNDEF 01429 END DO 01430 ENDIF 01431 ! 01432 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WG_METEO_FRANCE',1,ZHOOK_HANDLE) 01433 END SUBROUTINE READ_GRIB_WG_METEO_FRANCE 01434 !-------------------------------------------------------------------------------- 01435 ! ####################### 01436 SUBROUTINE READ_GRIB_WGI_METEO_FRANCE(HGRIB,KLUOUT,HINMODEL,PMASK,PFIELD,PD) 01437 ! ####################### 01438 ! 01439 USE MODD_GRID_GRIB, ONLY : NNI 01440 USE MODD_SURF_PAR, ONLY : XUNDEF 01441 ! 01442 IMPLICIT NONE 01443 ! 01444 !* dummy arguments 01445 ! --------------- 01446 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01447 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01448 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01449 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01450 REAL, DIMENSION(:,:), POINTER :: PFIELD ! field to initialize 01451 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 01452 ! 01453 !* local variables 01454 ! --------------- 01455 LOGICAL :: GISBA ! T: surface scheme in file is ISBA 01456 INTEGER(KIND=kindOfInt) :: IRET ! return code 01457 INTEGER :: ILTYPE ! type of level (Grib code table 3) 01458 INTEGER :: ILEV1 ! level definition 01459 INTEGER :: ILEV2 ! level definition 01460 REAL, DIMENSION(:), POINTER :: ZCLAY => NULL() ! clay fraction 01461 REAL, DIMENSION(:), POINTER :: ZSAND => NULL() ! sand fraction 01462 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() 01463 REAL, DIMENSION(:), ALLOCATABLE:: ZWSAT ! saturation 01464 INTEGER :: JL ! layer loop counter 01465 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01466 !------------------------------------------------------------------------------- 01467 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WGI_METEO_FRANCE',0,ZHOOK_HANDLE) 01468 ! 01469 ! 1. Search and read clay fraction if available 01470 ! ------------------------------------------ 01471 CALL READ_GRIB_SAND_CLAY_METEO_FRANCE(HGRIB,KLUOUT,HINMODEL,ZSAND,ZCLAY,GISBA) 01472 !------------------------------------------------------------------------------- 01473 IF (GISBA) THEN 01474 ALLOCATE(PFIELD(SIZE(ZSAND),2)) 01475 ALLOCATE(PD(SIZE(ZSAND),2)) 01476 ELSE 01477 ALLOCATE(PFIELD(NNI,2)) 01478 ALLOCATE(PD(NNI,2)) 01479 ENDIF 01480 ! 01481 PD(:,1) = 0. 01482 !------------------------------------------------------------------------------- 01483 ! 2. Read layer 1 soil ice 01484 ! --------------------- 01485 ILEV1 = 0 01486 IF (HINMODEL == 'ARPEGE' .OR. HINMODEL=='MOCAGE') THEN 01487 ILTYPE = 112 01488 ILEV2 = 1 01489 CALL READ_GRIB(HGRIB,KLUOUT,152,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 01490 ELSE 01491 ILTYPE = 105 01492 ILEV2 = 0 01493 CALL READ_GRIB(HGRIB,KLUOUT,139,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 01494 END IF 01495 ! 01496 IF (IRET == 0) THEN 01497 PFIELD(:,1) = ZFIELD(:) 01498 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_WGI_METEO_FRANCE: -> Soil ice level 1 is present' 01499 ELSE 01500 PFIELD(:,1) = 0. 01501 END IF 01502 !------------------------------------------------------------------------------- 01503 ! 3. Read layer 2 soil ice 01504 ! --------------------- 01505 IF (HINMODEL == 'ARPEGE' .OR. HINMODEL=='MOCAGE') THEN 01506 ILTYPE = 112 01507 ILEV1 = 0 01508 ILEV2 = 250 01509 CALL READ_GRIB(HGRIB,KLUOUT,152,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 01510 ELSE 01511 ILTYPE = 111 01512 ILEV1 = -1 01513 ILEV2 = -1 01514 CALL READ_GRIB(HGRIB,KLUOUT,139,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 01515 END IF 01516 ! 01517 IF (IRET == 0) THEN 01518 PFIELD(:,2) = ZFIELD(:) 01519 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_WGI_METEO_FRANCE: -> Soil ice level 2 is present' 01520 ELSE 01521 PFIELD(:,2) = 0. 01522 END IF 01523 !------------------------------------------------------------------------------- 01524 ! 4. Read layer 2 depth (ISBA only) 01525 ! ------------------------------ 01526 IF (GISBA) THEN 01527 IF (HINMODEL == 'ARPEGE' .OR. HINMODEL=='MOCAGE') THEN 01528 CALL READ_GRIB(HGRIB,KLUOUT,173,IRET,ZFIELD) 01529 ELSE 01530 CALL READ_GRIB(HGRIB,KLUOUT,130,IRET,ZFIELD) 01531 ENDIF 01532 IF (IRET /= 0) THEN 01533 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL ICE LEVEL 2 MISSING (READ_GRIB_WGI_METEO_FRANCE)') 01534 END IF 01535 PD(:,2) = ZFIELD(:) 01536 DEALLOCATE(ZFIELD) 01537 ELSE 01538 PD(:,2) = 2. 01539 END IF 01540 !------------------------------------------------------------------------------- 01541 ! 5. Compute relative humidity from units kg/m^2 01542 ! ------------------------------------------- 01543 IF (GISBA) THEN 01544 ! 01545 !* updates Wgi in m3/m3 01546 PFIELD(:,1) = PFIELD(:,1) / 10. 01547 PFIELD(:,2) = PFIELD(:,2) /(1000. * PD(:,2)) 01548 ! 01549 ALLOCATE (ZWSAT (NNI)) 01550 ZWSAT (:) = (-1.08*100. * ZSAND(:) + 494.305) * 0.001 01551 PFIELD(:,1) = PFIELD(:,1) / ZWSAT(:) 01552 PFIELD(:,2) = PFIELD(:,2) / ZWSAT(:) 01553 DEALLOCATE (ZWSAT) 01554 DEALLOCATE (ZSAND) 01555 DEALLOCATE (ZCLAY) 01556 ! 01557 ELSE ! Non ISBA 01558 ! 01559 PFIELD(:,1) = 0. 01560 PFIELD(:,2) = 0. 01561 ! 01562 END IF 01563 !-------------------------------------------------------------------------------- 01564 ! 6. Apply land mask 01565 ! --------------- 01566 IF (SIZE(PMASK)==SIZE(PFIELD,1)) THEN 01567 DO JL=1,SIZE(PFIELD,2) 01568 WHERE (PMASK(:)/=1.) PFIELD(:,JL) = XUNDEF 01569 END DO 01570 ENDIF 01571 ! 01572 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WGI_METEO_FRANCE',1,ZHOOK_HANDLE) 01573 END SUBROUTINE READ_GRIB_WGI_METEO_FRANCE 01574 !------------------------------------------------------------------- 01575 !------------------------------------------------------------------- 01576 ! ####################### 01577 SUBROUTINE READ_GRIB_TG_HIRLAM(HGRIB,KLUOUT,HINMODEL,PMASK,PTG,PDT) 01578 ! ####################### 01579 ! 01580 USE MODD_GRID_GRIB, ONLY : NNI 01581 USE MODD_SURF_PAR, ONLY : XUNDEF 01582 ! 01583 IMPLICIT NONE 01584 ! 01585 !* dummy arguments 01586 ! --------------- 01587 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01588 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01589 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01590 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01591 REAL, DIMENSION(:,:), POINTER :: PTG ! field to initialize 01592 REAL, DIMENSION(:,:), POINTER :: PDT ! thickness of each layer 01593 ! 01594 !* local variables 01595 ! --------------- 01596 INTEGER(KIND=kindOfInt) :: IRET ! return code 01597 INTEGER :: ILEV1 ! level definition 01598 INTEGER :: ILEV2 ! level definition 01599 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() 01600 INTEGER :: JL ! layer loop counter 01601 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01602 !-------------------------------------------------------------------------------- 01603 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TG_HIRLAM',0,ZHOOK_HANDLE) 01604 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_TG_HIRLAM: | Reading soil temperature' 01605 !-------------------------------------------------------------------------------- 01606 ! 1. Search and read level 1 (and its depth) 01607 ! ----------------------- 01608 ILEV1 = 904 01609 ILEV2 = -1 01610 CALL READ_GRIB(HGRIB,KLUOUT,11,IRET,ZFIELD,KLEV1=ILEV1,KLEV2=ILEV2) 01611 IF (IRET /= 0 ) THEN 01612 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL TEMPERATURE LEVEL 1 MISSING (READ_GRIB_TG_HIRLAM)') 01613 END IF 01614 ! 01615 ALLOCATE(PTG(SIZE(ZFIELD),3)) 01616 ALLOCATE(PDT(SIZE(ZFIELD),3)) 01617 PTG(:,1)= ZFIELD(:) 01618 PDT(:,1) = 0. 01619 !-------------------------------------------------------------------------------- 01620 ! 2. Deep soil temperature 01621 ! --------------------- 01622 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_TG_HIRLAM: | Reading deep soil temperature' 01623 ! 01624 ILEV1 = 954 01625 ILEV2 = -1 01626 CALL READ_GRIB(HGRIB,KLUOUT,11,IRET,ZFIELD,KLEV1=ILEV1,KLEV2=ILEV2) 01627 IF (IRET /= 0) THEN 01628 CALL ABOR1_SFX('MODE_READ_GRIB: DEEP SOIL TEMPERATURE MISSING (READ_GRIB_TG_HIRLAM)') 01629 END IF 01630 ! 01631 PTG(:,2)= ZFIELD(:) 01632 DEALLOCATE(ZFIELD) 01633 PDT(:,2) = 0.2 ! deep temperature depth assumed equal to 0.20m 01634 !-------------------------------------------------------------------------------- 01635 ! 4. Assumes uniform temperature profile below 01636 ! ----------------------------------------- 01637 PTG(:,3) = PTG(:,2) 01638 PDT(:,3) = 3. ! temperature profile down to 3m 01639 !-------------------------------------------------------------------------------- 01640 ! 5. Apply land mask 01641 ! --------------- 01642 IF (SIZE(PMASK)==SIZE(PTG,1)) THEN 01643 DO JL=1,SIZE(PTG,2) 01644 WHERE (PMASK(:)/=1.) PTG(:,JL) = XUNDEF 01645 END DO 01646 ENDIF 01647 ! 01648 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TG_HIRLAM',1,ZHOOK_HANDLE) 01649 END SUBROUTINE READ_GRIB_TG_HIRLAM 01650 !------------------------------------------------------------------- 01651 ! ####################### 01652 SUBROUTINE READ_GRIB_WG_HIRLAM(HGRIB,KLUOUT,HINMODEL,PMASK,PFIELD,PD) 01653 ! ####################### 01654 ! 01655 USE MODD_GRID_GRIB, ONLY : NNI 01656 USE MODD_SURF_PAR, ONLY : XUNDEF 01657 ! 01658 IMPLICIT NONE 01659 ! 01660 !* dummy arguments 01661 ! --------------- 01662 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01663 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01664 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01665 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01666 REAL, DIMENSION(:,:), POINTER :: PFIELD ! field to initialize 01667 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 01668 ! 01669 !* local variables 01670 ! --------------- 01671 INTEGER(KIND=kindOfInt) :: IRET ! return code 01672 INTEGER :: ILTYPE ! type of level (Grib code table 3) 01673 INTEGER :: ILEV1 ! level definition 01674 INTEGER :: ILEV2 ! level definition 01675 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() 01676 REAL, DIMENSION(:,:), ALLOCATABLE :: ZWG ! first water reservoir 01677 REAL, DIMENSION(:), ALLOCATABLE :: ZD ! Height of each layer 01678 INTEGER :: INLAYERDEEP! number of deep moisture layers 01679 REAL :: ZWWILT ! ECMWF wilting point 01680 REAL :: ZWFC ! ECMWF field capacity 01681 REAL :: ZWSAT ! ECMWF saturation 01682 INTEGER :: JL ! loop counter on layers 01683 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01684 !-------------------------------------------------------------------------------- 01685 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WG_HIRLAM',0,ZHOOK_HANDLE) 01686 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_WG_HIRLAM: | Reading soil moisture' 01687 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_WG_HIRLAM: | WARNING READING LOW VEGETATION TILE (NR 4) ONLY' 01688 ! 01689 ALLOCATE(ZD(2)) 01690 ! 01691 ! 1. Search and read level 1 (and its depth) 01692 ! ----------------------- 01693 ILTYPE=105 01694 ILEV1=904 01695 ILEV2=-1 01696 CALL READ_GRIB(HGRIB,KLUOUT,86,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 01697 IF (IRET /= 0 ) THEN 01698 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL MOISTURE LEVEL 1 MISSING (READ_GRIB_WG_HIRLAM)') 01699 END IF 01700 ! 01701 ALLOCATE(ZWG(SIZE(ZFIELD),2)) 01702 ZWG(:,1)=ZFIELD 01703 ! 01704 ZD(1) = 0.01 01705 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_WG_HIRLAM: | Level 1 height set to 0.01 m ' 01706 ! 01707 ZWG(:,1) = ZWG(:,1) / ZD(1) ! convert units to m3/m3 01708 !-------------------------------------------------------------------------------- 01709 ! 2. Search and read level 2 (and its depth) This level is optionnal 01710 ! ----------------------- 01711 ILTYPE=105 01712 ILEV1=954 01713 ILEV2=-1 01714 CALL READ_GRIB(HGRIB,KLUOUT,86,IRET,ZFIELD,KLTYPE=ILTYPE,KLEV1=ILEV1,KLEV2=ILEV2) 01715 IF (IRET /= 0 ) THEN 01716 CALL ABOR1_SFX('MODE_READ_GRIB: SOIL MOISTURE LEVEL 2 MISSING (READ_GRIB_WG_HIRLAM)') 01717 END IF 01718 ! 01719 ZWG(:,2)=ZFIELD ! already units m3/m3 01720 DEALLOCATE(ZFIELD) 01721 ! 01722 INLAYERDEEP = 2 01723 ZD(2) = 0.42 01724 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_WG_HIRLAM: | Level 2 height set to 0.42 m ' 01725 ! 01726 WRITE (KLUOUT,'(A)') 'WARNING MODE_READ_GRIB: ZWG3 AND ZWG4 SET TO 0. (READ_GRIB_WG_HIRLAM)' 01727 !-------------------------------------------------------------------------------- 01728 ! 3. Set water content profile and layer thicknesses 01729 ! ----------------------------------------------- 01730 CALL FILL_PFIELD(KLUOUT,'READ_GRIB_WG_HIRLAM',INLAYERDEEP,ZD,ZWG,PMASK,PFIELD,PD) 01731 DEALLOCATE(ZD) 01732 DEALLOCATE(ZWG) 01733 !-------------------------------------------------------------------------------- 01734 ! 4. Convert from specific humidity to relative humidity 01735 ! --------------------------------------------------- 01736 ! Compute model's constants 01737 ZWFC = 0.171 01738 ZWWILT = 0.086 01739 ! 01740 ! Then perform conversion 01741 DO JL=1,INLAYERDEEP 01742 WHERE (PFIELD(:,JL).NE.XUNDEF) PFIELD(:,JL) = (PFIELD(:,JL) - ZWWILT) / (ZWFC - ZWWILT) 01743 ENDDO 01744 ! 01745 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WG_HIRLAM',1,ZHOOK_HANDLE) 01746 END SUBROUTINE READ_GRIB_WG_HIRLAM 01747 !------------------------------------------------------------------- 01748 ! ####################### 01749 SUBROUTINE READ_GRIB_WGI_HIRLAM(HGRIB,KLUOUT,PFIELD,PD) 01750 ! ####################### 01751 ! 01752 USE MODD_GRID_GRIB, ONLY : NNI 01753 ! 01754 IMPLICIT NONE 01755 ! 01756 !* dummy arguments 01757 ! --------------- 01758 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01759 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01760 REAL, DIMENSION(:,:), POINTER :: PFIELD ! field to initialize 01761 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 01762 ! 01763 !* local variables 01764 ! --------------- 01765 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01766 !-------------------------------------------------------------------------------- 01767 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WGI_HIRLAM',0,ZHOOK_HANDLE) 01768 ! 01769 ALLOCATE (PFIELD(NNI,2)) 01770 ALLOCATE (PD (NNI,2)) 01771 PFIELD(:,:) = 0. 01772 ! 01773 PD (:,1) = 0. 01774 PD (:,2) = 1. 01775 ! 01776 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_WGI_HIRLAM',1,ZHOOK_HANDLE) 01777 END SUBROUTINE READ_GRIB_WGI_HIRLAM 01778 !------------------------------------------------------------------- 01779 !------------------------------------------------------------------- 01780 ! ####################### 01781 SUBROUTINE READ_GRIB_SNOW_VEG_AND_DEPTH(HGRIB,KLUOUT,HINMODEL,PMASK,PSNV,PSNVD) 01782 ! ####################### 01783 ! 01784 USE MODD_GRID_GRIB, ONLY : NNI 01785 USE MODD_SURF_PAR, ONLY : XUNDEF 01786 USE MODD_SNOW_PAR, ONLY : XRHOSMAX 01787 ! 01788 IMPLICIT NONE 01789 ! 01790 !* dummy arguments 01791 ! --------------- 01792 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01793 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01794 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01795 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01796 REAL, DIMENSION(:), OPTIONAL, POINTER :: PSNV ! field to initialize 01797 REAL, DIMENSION(:), OPTIONAL, POINTER :: PSNVD ! field to initialize 01798 ! 01799 !* local variables 01800 ! --------------- 01801 INTEGER(KIND=kindOfInt) :: IRET ! return code 01802 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() ! field to initialize 01803 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01804 !-------------------------------------------------------------------------------- 01805 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_SNOW_VEG_AND_DEPTH',0,ZHOOK_HANDLE) 01806 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_SNOW_VEG_AND_DEPTH: | Reading snow depth' 01807 ! 01808 SELECT CASE(HINMODEL) 01809 CASE('ECMWF ') 01810 CALL READ_GRIB(HGRIB,KLUOUT,141,IRET,ZFIELD) 01811 CASE('ARPEGE','ALADIN','MOCAGE','HIRLAM') 01812 CALL READ_GRIB(HGRIB,KLUOUT,66,IRET,ZFIELD) 01813 CASE DEFAULT 01814 CALL ABOR1_SFX('MODE_READ_GRIB:READ_GRIB_SNOW_VEG_AND_DEPTH: OPTION NOT SUPPORTED '//HINMODEL) 01815 END SELECT 01816 ! 01817 IF (IRET /= 0 ) THEN 01818 CALL ABOR1_SFX('MODE_READ_GRIB: SNOW AND VEG DEPTH MISSING (READ_GRIB_SNOW_VEG_AND_DEPTH)') 01819 END IF 01820 ! 01821 IF (PRESENT(PSNV)) THEN 01822 ALLOCATE(PSNV(SIZE(ZFIELD))) 01823 PSNV(:)=ZFIELD(:) 01824 IF (HINMODEL=='ECMWF ') PSNV(:) = PSNV(:) * XRHOSMAX 01825 IF (SIZE(PMASK)==SIZE(PSNV)) & 01826 WHERE (PMASK(:)/=1.) PSNV(:) = XUNDEF 01827 ENDIF 01828 ! 01829 IF (PRESENT(PSNVD)) THEN 01830 ALLOCATE(PSNVD(SIZE(ZFIELD))) 01831 PSNVD(:)=ZFIELD(:) 01832 IF (HINMODEL/='ECMWF ') PSNVD = PSNVD / XRHOSMAX 01833 IF (SIZE(PMASK)==SIZE(PSNVD)) & 01834 WHERE (PMASK(:)/=1.) PSNVD(:) = XUNDEF 01835 ENDIF 01836 ! 01837 DEALLOCATE(ZFIELD) 01838 ! 01839 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_SNOW_VEG_AND_DEPTH',1,ZHOOK_HANDLE) 01840 END SUBROUTINE READ_GRIB_SNOW_VEG_AND_DEPTH 01841 !------------------------------------------------------------------- 01842 !------------------------------------------------------------------- 01843 ! ####################### 01844 SUBROUTINE READ_GRIB_T_TEB(HGRIB,KLUOUT,HINMODEL,PTI,PMASK,PT,PD) 01845 ! ####################### 01846 ! 01847 USE MODD_GRID_GRIB, ONLY : NNI 01848 USE MODD_SURF_PAR, ONLY : XUNDEF 01849 ! 01850 IMPLICIT NONE 01851 ! 01852 !* dummy arguments 01853 ! --------------- 01854 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01855 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01856 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01857 REAL, INTENT(IN) :: PTI ! internal temperature 01858 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01859 REAL, DIMENSION(:,:), POINTER :: PT ! field to initialize 01860 REAL, DIMENSION(:,:), POINTER :: PD ! normalized grid 01861 ! 01862 !* local variables 01863 ! --------------- 01864 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() ! field to initialize 01865 INTEGER :: JL ! layer loop counter 01866 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01867 !-------------------------------------------------------------------------------- 01868 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_T_TEB',0,ZHOOK_HANDLE) 01869 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_T_TEB: | Reading temperature for buildings' 01870 ! 01871 CALL READ_GRIB_TS(HGRIB,KLUOUT,HINMODEL,PMASK,ZFIELD) 01872 ! 01873 ALLOCATE(PT(SIZE(ZFIELD),3)) 01874 ALLOCATE(PD(SIZE(ZFIELD),3)) 01875 ! 01876 PT(:,1) = ZFIELD 01877 DEALLOCATE(ZFIELD) 01878 PD(:,1) = 0. 01879 ! 01880 PT(:,2) = PTI 01881 PD(:,2) = 0.5 ! deep temperature depth assumed at half of wall or roof 01882 ! 01883 PT(:,3) = PTI 01884 PD(:,3) = 1. ! temperature at building interior 01885 ! 01886 IF (SIZE(PMASK)==SIZE(PT,1)) THEN 01887 DO JL=1,SIZE(PT,2) 01888 WHERE (PMASK(:)/=1.) PT(:,JL) = XUNDEF 01889 END DO 01890 ENDIF 01891 ! 01892 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_T_TEB',1,ZHOOK_HANDLE) 01893 END SUBROUTINE READ_GRIB_T_TEB 01894 !------------------------------------------------------------------- 01895 ! ####################### 01896 SUBROUTINE READ_GRIB_TF_TEB(HGRIB,KLUOUT,HINMODEL,PTI,PMASK,PTF,PD) 01897 ! ####################### 01898 ! 01899 USE MODD_SURF_PAR, ONLY : XUNDEF 01900 ! 01901 IMPLICIT NONE 01902 ! 01903 !* dummy arguments 01904 ! --------------- 01905 CHARACTER(LEN=*), INTENT(IN) :: HGRIB ! Grib file name 01906 INTEGER, INTENT(IN) :: KLUOUT ! logical unit of output listing 01907 CHARACTER(LEN=6), INTENT(IN) :: HINMODEL ! Grib originating model 01908 REAL, INTENT(IN) :: PTI ! internal temperature 01909 REAL, DIMENSION(:), INTENT(IN) :: PMASK ! grib land mask 01910 REAL, DIMENSION(:,:), POINTER :: PTF ! field to initialize 01911 REAL, DIMENSION(:,:), POINTER :: PD ! thickness of each layer 01912 ! 01913 !* local variables 01914 ! --------------- 01915 INTEGER(KIND=kindOfInt) :: IRET ! return code 01916 REAL, DIMENSION(:), POINTER :: ZFIELD => NULL() ! field to read 01917 INTEGER :: JL ! layer loop counter 01918 REAL(KIND=JPRB) :: ZHOOK_HANDLE 01919 !-------------------------------------------------------------------------------- 01920 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TF_TEB',0,ZHOOK_HANDLE) 01921 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_TF_TEB: | Reading temperature for building floor' 01922 ! 01923 ! 1. Deep soil temperature 01924 ! --------------------- 01925 ! 01926 WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_TF_TEB: | Reading deep soil temperature' 01927 ! 01928 CALL READ_GRIB_T2(HGRIB,KLUOUT,HINMODEL,PMASK,ZFIELD) 01929 ! 01930 ALLOCATE(PTF(SIZE(ZFIELD),3)) 01931 ALLOCATE(PD (SIZE(ZFIELD),3)) 01932 ! 01933 PTF(:,2) = ZFIELD(:) 01934 PD (:,2) = 0.5 ! deep temperature depth assumed at half of the floor 01935 ! 01936 DEALLOCATE(ZFIELD) 01937 ! 01938 ! 2. level 1 is the internal building temperature 01939 ! ----------------------- 01940 ! 01941 PTF(:,1) = PTI 01942 PD (:,1) = 0. 01943 ! 01944 ! 3. Assumes uniform temperature profile below 01945 ! ----------------------------------------- 01946 ! 01947 PTF(:,3) = PTF(:,2) 01948 PD (:,3) = 1. ! deep temperature value 01949 ! 01950 IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_TF_TEB',1,ZHOOK_HANDLE) 01951 !------------------------------------------------------------------------------ 01952 END SUBROUTINE READ_GRIB_TF_TEB 01953 !------------------------------------------------------------------- 01954 END MODULE MODE_READ_GRIB
1.8.0
