SURFEX v7.3: /home/dasprezs/EXPORT_v7_3/src/OFFLIN/oi

Go to the documentation of this file.
00001 SUBROUTINE OI_CONTROL (      &
00002  & LDINLINE,                 &
00003  & P__SURFTEMPERATURE,       &
00004  & P__SURFPREC_EAU_CON,      &
00005  & P__SURFPREC_EAU_GEC,      &
00006  & P__SURFPREC_NEI_CON,      &
00007  & P__SURFPREC_NEI_GEC,      &
00008  & P__ATMONEBUL_BASSE,       &
00009  & P__SURFXEVAPOTRANSP,      &
00010  & P__SURFFLU_LAT_MEVA,      &
00011  & P__SURFACCPLUIE,          &
00012  & P__SURFACCNEIGE,          &
00013  & P__SURFACCGRAUPEL,        &
00014  & P__CLSTEMPERATURE,        &
00015  & P__CLSHUMI_RELATIVE,      &
00016  & P__SURFIND_TERREMER,      &
00017  & P__SURFRESERV_NEIGE,      &
00018  & P__LON,                   &
00019  & P__LAT,                   &
00020  & LD_MASKEXT)
00021 
00022 ! ------------------------------------------------------------------------------------------
00023 !  *****************************************************************************************
00024 !
00025 !  Program to perform within SURFEX 
00026 !  a soil analysis for water content and temperature 
00027 !  using the Meteo-France optimum interpolation technique of Giard and Bazile (2000)
00028 !
00029 !  Derived from CANARI subroutines externalized by Lora Taseva (Dec. 2007)
00030 !
00031 !  Author : Jean-Francois Mahfouf (01/2008)
00032 !
00033 !  Modifications : 
00034 !   (05/2008)  : The I/O of this version follow the newly available LFI format in SURFEX  
00035 !   (01/2009)  : Read directly atmospheric FA files using XRD library instead of using "edf"
00036 !   (06/2009)  : Modifications to allow the assimilation of ASCAT superficial soil moisture
00037 !   (09/2010)  : More parameters to goto_surfex
00038 !   (03/2011)  : Initialization of ZEVAPTR (F.Bouyssel)
00039 !
00040 ! ******************************************************************************************
00041 ! ------------------------------------------------------------------------------------------
00042 USE MODD_TYPE_DATE_SURF
00043 USE MODD_CSTS,       ONLY : XDAY, XPI, XRHOLW, XLVTT, NDAYSEC
00044 USE MODD_SURF_PAR,   ONLY : XUNDEF
00045 USE MODD_ASSIM
00046 USE MODD_OL_FILEID
00047 
00048 USE MODD_SURFEX_OMP, ONLY : NINDX2_s=>NINDX2, NWORK, XWORK, XWORK2
00049 
00050 USE MODN_IO_OFFLINE, ONLY : CSURF_FILETYPE
00051 
00052 USE MODD_SURF_ATM_n, ONLY : CSEA,        CWATER,      CTOWN,      CNATURE,      &
00053                             XSEA,        XWATER,      XTOWN,      XNATURE,      &
00054                             NSIZE_SEA,   NSIZE_WATER, NSIZE_TOWN, NSIZE_NATURE, &
00055                             NR_SEA,      NR_WATER,    NR_TOWN,                  &
00056                             NR_NATURE,   XCOVER,      NDIM_FULL,  NSIZE_FULL,   &
00057                             NDIM_NATURE, NDIM_SEA,    NDIM_WATER, NDIM_TOWN 
00058 
00059 #ifdef LFI 
00060 USE MODD_IO_SURF_LFI,ONLY : CFILEIN_LFI, CFILEOUT_LFI, CFILEPGD_LFI, CFILEIN_LFI_SAVE
00061 #endif
00062 #ifdef FA
00063 USE MODD_IO_SURF_FA, ONLY : CFILEIN_FA, CFILEIN_FA_SAVE, CDNOMC, &
00064                             NDGUX,  NDLUX,  PERPK,  PELON0, PELAT0, &
00065                             PEDELX, PEDELY, PELON1, PELAT1, PEBETA
00066 #endif 
00067 #ifdef ARO 
00068 USE MODD_IO_SURF_ARO,ONLY : NGPTOT, NGPTOT_CAP, NPROMA, NINDX1, NINDX2, NBLOCK, NKPROMA, &
00069                             LWRITE, LCOUNTW, LFMWRIT, XGPGW, YSURFEX_CACHE_OUT,          &
00070                             SURFEX_FIELD_BUF_PREALLOC, SURFEX_FIELD_BUF_SET_RECORD,      &
00071                             NCOUNTW, NCOUNTW_TOT 
00072 #endif
00073 
00074 USE MODI_READ_ALL_NAMELISTS
00075 USE MODI_GOTO_SURFEX
00076 USE MODI_ALLOC_SURFEX
00077 USE MODI_INI_DATA_COVER
00078 USE MODI_INIT_IO_SURF_n
00079 USE MODI_READ_SURF
00080 USE MODI_TRANS_CHAINE
00081 USE MODI_SET_SURFEX_FILEIN
00082 USE MODI_GET_SIZE_FULL_n
00083 USE MODI_READ_COVER_n
00084 USE MODI_CONVERT_COVER_FRAC
00085 USE MODI_GET_1D_MASK
00086 USE MODI_END_IO_SURF_n
00087 USE MODI_IO_BUFF_CLEAN_n
00088 USE MODI_OI_BC_SOIL_MOISTURE
00089 USE MODI_OI_LATLON_CONF_PROJ
00090 USE MODI_OI_CACSTS
00091 USE MODI_OI_HOR_EXTRAPOL_SURF
00092 USE MODI_FLAG_UPDATE
00093 USE MODI_WRITE_SURF
00094 
00095 USE YOMHOOK ,ONLY : LHOOK, DR_HOOK
00096 USE PARKIND1  ,ONLY : JPRB
00097 
00098 IMPLICIT NONE
00099 
00100 LOGICAL, INTENT (IN) :: LDINLINE
00101 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFTEMPERATURE
00102 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFPREC_EAU_CON
00103 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFPREC_EAU_GEC
00104 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFPREC_NEI_CON
00105 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFPREC_NEI_GEC
00106 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__ATMONEBUL_BASSE
00107 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFXEVAPOTRANSP
00108 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFFLU_LAT_MEVA
00109 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFACCPLUIE
00110 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFACCNEIGE
00111 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFACCGRAUPEL
00112 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__CLSTEMPERATURE
00113 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__CLSHUMI_RELATIVE
00114 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFIND_TERREMER
00115 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__SURFRESERV_NEIGE
00116 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__LON
00117 REAL(KIND=JPRB), OPTIONAL, DIMENSION (:) ::  P__LAT
00118 LOGICAL,         OPTIONAL, DIMENSION (:) ::  LD_MASKEXT
00119 
00120 INTEGER :: IGPCOMP
00121 INTEGER :: IDAT
00122 
00123  CHARACTER(LEN=28) :: YNAMELIST = 'OPTIONS.nam                 '
00124 
00125 !    Declarations of local variables
00126 
00127  CHARACTER(LEN=6)            :: YPROGRAM
00128  CHARACTER(LEN=6), PARAMETER :: YPROGRAM2 = 'FA    '
00129  CHARACTER(LEN=2)            :: CMONTH
00130 INTEGER                     :: IYEAR                ! current year (UTC)
00131 INTEGER                     :: IMONTH               ! current month (UTC)
00132 INTEGER                     :: IDAY                 ! current day (UTC)
00133 INTEGER                     :: NSSSSS               ! current time since start of the run (s)
00134 INTEGER                     :: IRESP                ! return code
00135 TYPE (DATE_TIME)            :: TTIME                ! Current date and time
00136 INTEGER                     :: ISIZE
00137 INTEGER                     :: ISIZE1
00138 LOGICAL                     :: LLKEEPEXTZONE
00139 
00140 ! Arrays for soil OI analysis
00141 REAL, DIMENSION (:,:), ALLOCATABLE :: PWS, PWP, PTS, PTP, PTL, PSNS, PRSMIN, PD2, PLAI, PVEG
00142 REAL, DIMENSION (:),   ALLOCATABLE :: PSST, PSAB, PARG, PLAT, PLON, PTCLS, PHCLS, PUCLS, PVCLS,  
00143                                       PEVAP, PEVAPTR, PT2M_O, PHU2M_O, PTS_O, ZT2INC, ZH2INC,   &
00144                                       ZWS, ZWP, ZTL, ZTS, ZTP, ZTCLS, ZHCLS, ZUCLS, ZVCLS,      &
00145                                       PSSTC, PWPINC1, PWPINC2, PWPINC3, PT2MBIAS, PH2MBIAS,     &
00146                                       PRRCN, PRRCL, PRRSN, PRRSL, PATMNEB, PITM, PALBF, PEMISF, &
00147                                       PZ0F, PIVEG, PZ0H, PTSC, PTPC, PWSC, PWPC, PSNC, ZEVAP,   &
00148                                       ZEVAPTR, PGELAT, PGELAM, PGEMU, ZWSINC, ZWPINC, ZTSINC,   &
00149                                       ZTPINC, ZTLINC, ZSNINC, ZSNS, ZPX, ZPY, PSM_O, PSIG_SMO,  &
00150                                       PLSM_O, PWS_O, ZWGINC, PLST, PTRD3, ZSST, ZLST, ZALT
00151 REAL, DIMENSION (:),   ALLOCATABLE :: ZSST1, ZLST1, PSST1, PLST1, PLAT1, PLON1, ZALT1
00152 
00153 INTEGER                            :: IVERSION, IBUGFIX
00154 INTEGER                            :: J,J1
00155  CHARACTER(LEN=10)                  :: YVAR    ! Name of the prognostic variable (in LFI file)
00156  CHARACTER(LEN=100)                 :: YPREFIX ! Prefix of the prognostic variable  (in LFI file)
00157 INTEGER                            :: ILUOUT  ! ascii output unit number
00158 INTEGER                            :: INOBS   ! number of observations
00159 
00160 REAL                               :: PLAT0,PLON0,PRPK,PLATOR,PLONOR,DELX,DELY,PBETA,ZTHRES
00161 
00162 LOGICAL, DIMENSION(:), ALLOCATABLE :: OINTERP_LST, OINTERP_SST
00163 LOGICAL, DIMENSION(:), ALLOCATABLE :: OINTERP_LST1, OINTERP_SST1
00164 
00165 REAL(KIND=JPRB) :: ZHOOK_HANDLE
00166 ! ----------------------------------------------------------------------------------
00167 IF (LHOOK) CALL DR_HOOK ('OI_CONTROL', 0, ZHOOK_HANDLE)
00168 
00169 PRINT *,'--------------------------------------------------------------------------'
00170 PRINT *,'|                                                                        |'
00171 PRINT *,'|                             ENTER OI_ASSIM                             |'
00172 PRINT *,'|                                                                        |'
00173 PRINT *,'--------------------------------------------------------------------------'
00174 
00175  CALL READ_ALL_NAMELISTS(CSURF_FILETYPE,'ALL',.FALSE.)
00176 
00177 IF (LDINLINE) THEN
00178 
00179   YPROGRAM = 'AROME'
00180 #ifdef ARO   
00181   IGPCOMP = MIN (NGPTOT, NGPTOT_CAP)
00182   
00183   NBLOCK   = 1
00184   NINDX1   = 1 + (NBLOCK - 1) * NPROMA
00185   NINDX2   = MIN (NBLOCK * NPROMA, IGPCOMP)
00186   NKPROMA  = NINDX2 - NINDX1 + 1
00187   CALL GOTO_SURFEX (NBLOCK, .TRUE.)
00188 #endif
00189 
00190 ELSE
00191 
00192   YPROGRAM = 'LFI'
00193   CALL ALLOC_SURFEX(1)
00194   CALL GOTO_SURFEX(1,.TRUE.)
00195 
00196 ENDIF
00197 
00198 ILUOUT = 111
00199 LLKEEPEXTZONE = .FALSE.
00200 
00201 !   Update some constants dependant from NACVEG
00202 
00203 !  scaling of soil moisture increments when assimilation window is different from 6 hours
00204 RSCALDW = REAL(NECHGU)/6.0
00205 !  half assimilation window in sec
00206 ITRAD   = NECHGU*1800
00207 
00208  CALL INI_DATA_COVER
00209 
00210 !   File handling definition
00211 
00212 IF (.NOT. LDINLINE) THEN
00213 #ifdef LFI
00214   CFILEPGD_LFI = 'PGD'
00215   CFILEIN_LFI = 'PREP'        ! input PREP file (surface fields) 
00216   CFILEIN_LFI_SAVE = CFILEIN_LFI
00217 #endif
00218 ENDIF
00219 
00220 !   Read grid dimension for allocation
00221 
00222  CALL INIT_IO_SURF_n(YPROGRAM,'FULL  ','SURF  ','READ ')
00223 
00224 !   Find current time
00225 
00226  CALL READ_SURF(YPROGRAM,'DTCUR',TTIME,IRESP)
00227 
00228 !   Time initializations
00229 
00230 IYEAR  = TTIME%TDATE%YEAR
00231 IMONTH = TTIME%TDATE%MONTH
00232 IDAY   = TTIME%TDATE%DAY
00233 NSSSSS = TTIME%TIME
00234 IF (NSSSSS > NDAYSEC) NSSSSS = NSSSSS - NDAYSEC
00235  CALL TRANS_CHAINE(CMONTH,IMONTH,2)
00236 
00237 !   Reading grid characteristics to perform nature mask
00238 
00239  CALL END_IO_SURF_n(YPROGRAM)
00240  CALL SET_SURFEX_FILEIN(YPROGRAM,'PGD ') ! change input file name to pgd name
00241  CALL INIT_IO_SURF_n(YPROGRAM,'FULL  ','SURF  ','READ ')
00242 
00243  CALL READ_SURF(YPROGRAM,'SEA   ',CSEA   ,IRESP)
00244  CALL READ_SURF(YPROGRAM,'WATER ',CWATER ,IRESP)
00245  CALL READ_SURF(YPROGRAM,'NATURE',CNATURE,IRESP)
00246  CALL READ_SURF(YPROGRAM,'TOWN  ',CTOWN  ,IRESP)
00247 
00248  CALL READ_SURF(YPROGRAM,'DIM_FULL  ',NDIM_FULL,  IRESP)
00249  CALL READ_SURF(YPROGRAM,'DIM_SEA   ',NDIM_SEA,   IRESP)
00250  CALL READ_SURF(YPROGRAM,'DIM_NATURE',NDIM_NATURE,IRESP)
00251  CALL READ_SURF(YPROGRAM,'DIM_WATER ',NDIM_WATER, IRESP)
00252  CALL READ_SURF(YPROGRAM,'DIM_TOWN  ',NDIM_TOWN,  IRESP)
00253 
00254 NINDX2_s = NDIM_FULL
00255 ALLOCATE(NWORK(NDIM_FULL))
00256 ALLOCATE(XWORK(NDIM_FULL))
00257 ALLOCATE(XWORK2(NDIM_FULL,10))
00258 
00259 !   Get total dimension of domain (excluding extension zone)
00260 
00261  CALL GET_SIZE_FULL_n(YPROGRAM,NDIM_FULL,NSIZE_FULL)
00262 
00263 IF (LDINLINE) THEN
00264   ISIZE = NSIZE_FULL
00265 ELSE
00266   ISIZE = NDIM_FULL
00267 ENDIF
00268 
00269 ALLOCATE (PSAB(ISIZE)) 
00270 ALLOCATE (PARG(ISIZE))
00271 ALLOCATE (ZALT(ISIZE))
00272 
00273  CALL READ_SURF(YPROGRAM,'SAND',      PSAB,  IRESP)
00274  CALL READ_SURF(YPROGRAM,'CLAY',      PARG,  IRESP)
00275  CALL READ_SURF(YPROGRAM,'ZS',        ZALT,  IRESP)
00276 
00277  CALL READ_COVER_n(YPROGRAM)
00278 
00279 !   Perform masks (only nature used)
00280 
00281 ALLOCATE(XSEA   (ISIZE))
00282 ALLOCATE(XNATURE(ISIZE))
00283 ALLOCATE(XWATER (ISIZE))
00284 ALLOCATE(XTOWN  (ISIZE))
00285 
00286  CALL CONVERT_COVER_FRAC(XCOVER,XSEA,XNATURE,XTOWN,XWATER)
00287 
00288 NSIZE_NATURE = COUNT(XNATURE(:) > 0.0)
00289 NSIZE_TOWN   = COUNT(XTOWN(:)   > 0.0)
00290 NSIZE_WATER  = COUNT(XWATER(:)  > 0.0)
00291 NSIZE_SEA    = COUNT(XSEA(:)    > 0.0)
00292 
00293 ALLOCATE(NR_NATURE (NSIZE_NATURE))
00294 ALLOCATE(NR_TOWN   (NSIZE_TOWN  ))
00295 ALLOCATE(NR_WATER  (NSIZE_WATER ))
00296 ALLOCATE(NR_SEA    (NSIZE_SEA   ))
00297 
00298  CALL GET_1D_MASK( NSIZE_SEA,    ISIZE, XSEA   , NR_SEA   )
00299  CALL GET_1D_MASK( NSIZE_WATER,  ISIZE, XWATER , NR_WATER )
00300  CALL GET_1D_MASK( NSIZE_TOWN,   ISIZE, XTOWN  , NR_TOWN  )
00301  CALL GET_1D_MASK( NSIZE_NATURE, ISIZE, XNATURE, NR_NATURE)
00302 
00303 ! Allocate arrays
00304 
00305 ALLOCATE (PWS(ISIZE,1))
00306 ALLOCATE (PWP(ISIZE,1))
00307 ALLOCATE (PTS(ISIZE,1))
00308 ALLOCATE (PTP(ISIZE,1))
00309 ALLOCATE (PTL(ISIZE,1))
00310 ALLOCATE (PSST(ISIZE))
00311 ALLOCATE (PSNS(ISIZE,1))
00312 ALLOCATE (PLAI(ISIZE,1))
00313 ALLOCATE (PVEG(ISIZE,1))
00314 ALLOCATE (PRSMIN(ISIZE,1))
00315 ALLOCATE (PD2(ISIZE,1))
00316 ALLOCATE (PTCLS(ISIZE))
00317 ALLOCATE (PHCLS(ISIZE))
00318 ALLOCATE (PUCLS(ISIZE))
00319 ALLOCATE (PVCLS(ISIZE))
00320 ALLOCATE (PEVAP(ISIZE))
00321 ALLOCATE (PLST(ISIZE))
00322 ALLOCATE (PTRD3(ISIZE))
00323 
00324 ALLOCATE (OINTERP_LST(ISIZE))
00325 ALLOCATE (OINTERP_SST(ISIZE))
00326 ALLOCATE (ZLST(ISIZE))
00327 ALLOCATE (ZSST(ISIZE))
00328 
00329 !  Read prognostic variables
00330 
00331  CALL END_IO_SURF_n(YPROGRAM)
00332  CALL SET_SURFEX_FILEIN(YPROGRAM,'PREP') ! change input file name to pgd name
00333  CALL INIT_IO_SURF_n(YPROGRAM,'FULL  ','SURF  ','READ ')
00334 
00335 IF (NSIZE_NATURE>0 .AND. CNATURE/='NONE') THEN
00336   CALL READ_SURF(YPROGRAM,'WG1',       PWS,   IRESP)
00337   CALL READ_SURF(YPROGRAM,'WG2',       PWP,   IRESP)
00338   CALL READ_SURF(YPROGRAM,'TG1',       PTS,   IRESP)
00339   CALL READ_SURF(YPROGRAM,'TG2',       PTP,   IRESP)
00340   CALL READ_SURF(YPROGRAM,'WGI2',      PTL,   IRESP)
00341 
00342   CALL READ_SURF(YPROGRAM,'VERSION',IVERSION,IRESP)
00343   CALL READ_SURF(YPROGRAM,'BUG',IBUGFIX,IRESP)
00344   IF (IVERSION>7 .OR. IVERSION==7 .AND. IBUGFIX>=3) THEN
00345     CALL READ_SURF(YPROGRAM,'WSN_VEG1',PSNS,  IRESP)
00346   ELSE
00347     CALL READ_SURF(YPROGRAM,'WSNOW_VEG1',PSNS,  IRESP)
00348   ENDIF
00349 ENDIF
00350 
00351 IF (NSIZE_SEA>0 .AND. CSEA/='NONE') THEN
00352   CALL READ_SURF(YPROGRAM,'SST',       PSST,  IRESP)
00353 ENDIF
00354 
00355 IF (NSIZE_WATER>0 .AND. CWATER/='NONE') THEN
00356   CALL READ_SURF(YPROGRAM,'TS_WATER',  PLST,  IRESP)
00357 ENDIF
00358 
00359 IF (NSIZE_TOWN>0 .AND. CTOWN/='NONE' .AND. LAROME) THEN
00360   CALL READ_SURF(YPROGRAM,'VERSION',IVERSION,IRESP)
00361   CALL READ_SURF(YPROGRAM,'BUG',IBUGFIX,IRESP)
00362   IF (IVERSION>7 .OR. IVERSION==7 .AND. IBUGFIX>=3) THEN
00363     CALL READ_SURF(YPROGRAM,'TROAD3',   PTRD3,  IRESP)
00364   ELSE
00365     CALL READ_SURF(YPROGRAM,'T_ROAD3',   PTRD3,  IRESP)
00366   ENDIF
00367 ELSE
00368   PTRD3(:) = XUNDEF
00369 ENDIF
00370 
00371  CALL READ_SURF(YPROGRAM,'T2M',       PTCLS, IRESP)
00372  CALL READ_SURF(YPROGRAM,'HU2M',      PHCLS, IRESP)
00373  CALL READ_SURF(YPROGRAM,'ZON10M',    PUCLS, IRESP)
00374  CALL READ_SURF(YPROGRAM,'MER10M',    PVCLS, IRESP)
00375 
00376 ! Read constant surface fields
00377 
00378  CALL READ_SURF(YPROGRAM,'RSMIN',     PRSMIN,IRESP)
00379  CALL READ_SURF(YPROGRAM,'DG2',       PD2,   IRESP)
00380  CALL READ_SURF(YPROGRAM,'LAI',       PLAI,  IRESP)
00381  CALL READ_SURF(YPROGRAM,'VEG',       PVEG,  IRESP)
00382 
00383 IF (LPRINT) THEN
00384   J = NR_NATURE(1)
00385   PRINT *,'value in PREP file => WG1       ',PWS(J,1)
00386   PRINT *,'value in PREP file => WG2       ',PWP(J,1)
00387   PRINT *,'value in PREP file => TG1       ',PTS(J,1)
00388   PRINT *,'value in PREP file => TG2       ',PTP(J,1)
00389   PRINT *,'value in PREP file => WGI2      ',PTL(J,1)
00390   PRINT *,'value in PREP file => WSNOW_VEG1',PSNS(J,1)
00391   PRINT *,'value in PREP file => SST       ',PSST(J)
00392   PRINT *,'value in PREP file => LAI       ',PLAI(J,1)
00393   PRINT *,'value in PREP file => VEG       ',PVEG(J,1)
00394   PRINT *,'value in PREP file => RSMIN     ',PRSMIN(J,1)
00395   PRINT *,'value in PREP file => DATA_DG2  ',PD2(J,1)
00396   PRINT *,'value in PREP file => SAND      ',PSAB(J)
00397   PRINT *,'value in PREP file => CLAY      ',PARG(J)
00398   PRINT *,'value in PREP file => ZS        ',ZALT(J)
00399 ENDIF
00400 
00401  CALL END_IO_SURF_n(YPROGRAM)
00402  CALL IO_BUFF_CLEAN_n
00403 
00404 !  Interface (allocate arrays)
00405 
00406 ALLOCATE (PLAT(ISIZE))
00407 ALLOCATE (PLON(ISIZE))
00408 ALLOCATE (ZPX(ISIZE))
00409 ALLOCATE (ZPY(ISIZE))
00410 ALLOCATE (PEVAPTR(ISIZE))
00411 ALLOCATE (ZWP(ISIZE))
00412 ALLOCATE (ZWS(ISIZE))
00413 ALLOCATE (ZTL(ISIZE))
00414 ALLOCATE (ZTS(ISIZE))
00415 ALLOCATE (ZTP(ISIZE))
00416 ALLOCATE (ZSNS(ISIZE))
00417 ALLOCATE (ZTCLS(ISIZE))
00418 ALLOCATE (ZHCLS(ISIZE))
00419 ALLOCATE (ZUCLS(ISIZE))
00420 ALLOCATE (ZVCLS(ISIZE))
00421 ALLOCATE (PSSTC(ISIZE))
00422 ALLOCATE (PWPINC1(ISIZE))
00423 ALLOCATE (PWPINC2(ISIZE))
00424 ALLOCATE (PWPINC3(ISIZE))
00425 ALLOCATE (PT2MBIAS(ISIZE))
00426 ALLOCATE (PH2MBIAS(ISIZE))
00427 ALLOCATE (PRRCN(ISIZE))
00428 ALLOCATE (PRRCL(ISIZE))
00429 ALLOCATE (PRRSN(ISIZE))
00430 ALLOCATE (PRRSL(ISIZE))
00431 ALLOCATE (PATMNEB(ISIZE))
00432 ALLOCATE (PITM(ISIZE))
00433 ALLOCATE (PALBF(ISIZE))
00434 ALLOCATE (PEMISF(ISIZE))
00435 ALLOCATE (PZ0F(ISIZE))
00436 ALLOCATE (PIVEG(ISIZE))
00437 ALLOCATE (PZ0H(ISIZE))
00438 ALLOCATE (PTSC(ISIZE))
00439 ALLOCATE (PTPC(ISIZE))
00440 ALLOCATE (PWSC(ISIZE))
00441 ALLOCATE (PWPC(ISIZE))
00442 ALLOCATE (PSNC(ISIZE))
00443 ALLOCATE (ZEVAP(ISIZE))
00444 ALLOCATE (ZEVAPTR(ISIZE))
00445 ALLOCATE (PGELAT(ISIZE))
00446 ALLOCATE (PGELAM(ISIZE))
00447 ALLOCATE (PGEMU(ISIZE))
00448 ALLOCATE (PT2M_O(ISIZE))
00449 ALLOCATE (PHU2M_O(ISIZE))
00450 ALLOCATE (PTS_O(ISIZE))
00451 ALLOCATE (PSM_O(ISIZE))
00452 ALLOCATE (PSIG_SMO(ISIZE))
00453 ALLOCATE (PLSM_O(ISIZE))
00454 ALLOCATE (PWS_O(ISIZE))
00455 ALLOCATE (ZWGINC(ISIZE))
00456 
00457 IF (.NOT. LDINLINE) THEN
00458 
00459 !  Read atmospheric forecast fields from FA files 
00460 #ifdef FA
00461   CFILEIN_FA = 'FG_OI_MAIN'        ! input forecast
00462   CFILEIN_FA_SAVE  = CFILEIN_FA
00463 #endif
00464 !  Open FA file (LAM version with extension zone)
00465   CALL INIT_IO_SURF_n(YPROGRAM2,'EXTZON','SURF  ','READ ')
00466 ENDIF
00467 
00468 !  Read model forecast quantities
00469 
00470 IF (LAROME) THEN
00471   IF (LDINLINE) THEN
00472     PRRSL(:) = P__SURFACCPLUIE   (1:ISIZE)
00473     PRRSN(:) = P__SURFACCNEIGE   (1:ISIZE)
00474     PRRCN(:) = P__SURFACCGRAUPEL (1:ISIZE)
00475   ELSE
00476     CALL READ_SURF(YPROGRAM2,'SURFACCPLUIE',    PRRSL  ,IRESP)
00477     CALL READ_SURF(YPROGRAM2,'SURFACCNEIGE',    PRRSN  ,IRESP)
00478     CALL READ_SURF(YPROGRAM2,'SURFACCGRAUPEL',  PRRCN  ,IRESP)
00479   ENDIF
00480   PRRCL(:) = 0.0
00481 !   CALL READ_SURF(YPROGRAM2,'SURFIND.VEG.DOMI',PIVEG  ,IRESP) 
00482   PIVEG(:) = 0.0
00483 ELSE
00484   IF (LDINLINE) THEN
00485     PRRCL(:) = P__SURFPREC_EAU_CON (1:ISIZE)
00486     PRRSL(:) = P__SURFPREC_EAU_GEC (1:ISIZE)
00487     PRRCN(:) = P__SURFPREC_NEI_CON (1:ISIZE)
00488     PRRSN(:) = P__SURFPREC_NEI_GEC (1:ISIZE)
00489   ELSE
00490     CALL READ_SURF(YPROGRAM2,'SURFPREC.EAU.CON',PRRCL  ,IRESP)
00491     CALL READ_SURF(YPROGRAM2,'SURFPREC.EAU.GEC',PRRSL  ,IRESP)
00492     CALL READ_SURF(YPROGRAM2,'SURFPREC.NEI.CON',PRRCN  ,IRESP)
00493     CALL READ_SURF(YPROGRAM2,'SURFPREC.NEI.GEC',PRRSN  ,IRESP)
00494   ENDIF
00495   PIVEG(:) = 0.0
00496 ENDIF
00497 IF (LDINLINE) THEN
00498   PATMNEB(:) = P__ATMONEBUL_BASSE  (1:ISIZE)
00499   PITM(:)    = P__SURFIND_TERREMER (1:ISIZE)
00500   PEVAP(:)   = P__SURFFLU_LAT_MEVA (1:ISIZE)
00501 ELSE
00502   CALL READ_SURF(YPROGRAM2,'ATMONEBUL.BASSE ',PATMNEB,IRESP)
00503   CALL READ_SURF(YPROGRAM2,'SURFIND.TERREMER',PITM   ,IRESP) 
00504   CALL READ_SURF(YPROGRAM2,'SURFFLU.LAT.MEVA',PEVAP  ,IRESP) ! accumulated fluxes (not available in LFI)
00505 ENDIF
00506 IF (.NOT.LALADSURF) THEN
00507   IF (LDINLINE) THEN
00508     PEVAPTR(:) = P__SURFXEVAPOTRANSP (1:ISIZE)
00509   ELSE
00510     CALL READ_SURF(YPROGRAM2,'SURFXEVAPOTRANSP',PEVAPTR,IRESP) ! not in ALADIN SURFEX
00511   ENDIF
00512 ELSE
00513   PEVAPTR(:) = 0.0
00514 ENDIF
00515 
00516 IF (.NOT. LDINLINE) THEN
00517 !  Close FA file
00518   CALL END_IO_SURF_n(YPROGRAM2)
00519   CALL IO_BUFF_CLEAN_n
00520   PRINT *,'READ FG_OI_MAIN OK'
00521 ENDIF
00522 
00523 IF (.NOT. LDINLINE) THEN
00524 !  Define FA file name for CANARI analysis
00525 #ifdef FA
00526   CFILEIN_FA = 'CANARI'        ! input CANARI analysis
00527   CFILEIN_FA_SAVE  = CFILEIN_FA
00528 #endif
00529 !  Open FA file 
00530   CALL INIT_IO_SURF_n(YPROGRAM2,'EXTZON','SURF  ','READ ')
00531 ENDIF
00532 
00533 IF (LDINLINE) THEN
00534   PT2M_O(:)  = P__CLSTEMPERATURE   (1:ISIZE)
00535   PHU2M_O(:) = P__CLSHUMI_RELATIVE (1:ISIZE)
00536   PTS_O(:)   = P__SURFTEMPERATURE  (1:ISIZE)
00537 ELSE
00538 !  Read CANARI analysis
00539   CALL READ_SURF(YPROGRAM2,'CLSTEMPERATURE  ',PT2M_O ,IRESP)
00540   CALL READ_SURF(YPROGRAM2,'CLSHUMI.RELATIVE',PHU2M_O,IRESP)
00541   CALL READ_SURF(YPROGRAM2,'SURFTEMPERATURE ',PTS_O  ,IRESP)
00542 ENDIF
00543 
00544 IF (.NOT. LDINLINE) THEN
00545 !  Close CANARI file
00546   CALL END_IO_SURF_n(YPROGRAM2)
00547   CALL IO_BUFF_CLEAN_n
00548   PRINT *,'READ CANARI OK'
00549 ENDIF
00550 
00551 !  Read ASCAT SM observations (in percent)
00552 
00553 INOBS = 0
00554 IF (LOBSWG) THEN
00555   OPEN(UNIT=111,FILE='ASCAT_SM.DAT')
00556   DO J=1,NDIM_FULL
00557     READ(111,*,END=990) PSM_O(J),PSIG_SMO(J),PLSM_O(J)
00558     IF (PLSM_O(J) < 1.0)          PSM_O(J) = 999.0 ! data rejection if not on land
00559     IF (PSIG_SMO(J) > SIGWGO_MAX) PSM_O(J) = 999.0 ! data rejection of error too large
00560     IF (PSM_O(J) /= 999.0) INOBS = INOBS + 1
00561   ENDDO
00562 990 CONTINUE
00563   CLOSE(UNIT=111)
00564   PRINT *,'READ ASCAT SM OK'
00565 ELSE
00566   PSM_O(:)    = 999.0
00567   PSIG_SMO(:) = 999.0
00568   PLSM_O(:)   = 0.0
00569 ENDIF
00570 PRINT *,' NUMBER OF ASCAT OBSERVATIONS AFTER INITIAL CHECKS  :: ',INOBS
00571 INOBS = 0
00572 
00573 ! Perform bias correction of SM observations
00574 
00575  CALL OI_BC_SOIL_MOISTURE(ISIZE,PSM_O,PSAB,PWS_O)
00576 
00577 IF (.NOT. LDINLINE) THEN
00578 !  Define FA file name for surface climatology
00579 #ifdef FA
00580   CFILEIN_FA = 'clim_isba'               ! input climatology
00581   CFILEIN_FA_SAVE  = CFILEIN_FA
00582   CDNOMC     = 'climat'                  ! new frame name
00583 #endif
00584 !  Open FA file 
00585   CALL INIT_IO_SURF_n(YPROGRAM2,'EXTZON','SURF  ','READ ')
00586 ENDIF
00587 
00588 IF (LDINLINE) THEN
00589   PSNC(:) = P__SURFRESERV_NEIGE (1:ISIZE)
00590 ELSE
00591 !  Read climatology file (snow water equivalent)
00592   CALL READ_SURF(YPROGRAM2,'SURFRESERV.NEIGE',PSNC  ,IRESP)
00593 ENDIF
00594 
00595 IF (.NOT. LDINLINE) THEN
00596 !  Close climatology file
00597   CALL END_IO_SURF_n(YPROGRAM2)
00598   CALL IO_BUFF_CLEAN_n
00599   PRINT *,'READ CLIMATOLOGY OK'
00600 ENDIF
00601 
00602 IF (.NOT. LDINLINE) THEN
00603 #ifdef FA
00604   PLAT0  = PELAT0 
00605   PLON0  = PELON0 
00606   PLATOR = PELAT1 
00607   PLONOR = PELON1 
00608   PRPK   = PERPK  
00609   PBETA  = PEBETA 
00610   DELX   = PEDELX 
00611   DELY   = PEDELY 
00612   IF (PLONOR > 180.0) PLONOR = PLONOR - 360.0
00613   IF (PLON0  > 180.0) PLON0  = PLON0  - 360.0
00614   DO J=1,NDGUX
00615     DO J1=1,NDLUX
00616       ZPX((J-1)*NDLUX + J1) = DELX*REAL(J1-1)
00617       ZPY((J-1)*NDLUX + J1) = DELY*REAL(J-1)
00618     ENDDO
00619   ENDDO
00620 #endif
00621   CALL OI_LATLON_CONF_PROJ(ISIZE,PLAT0,PLON0,PRPK,PBETA,PLATOR,PLONOR,ZPX,ZPY,PLAT,PLON)
00622 ELSE
00623   PLAT(:) = P__LAT (1:ISIZE)
00624   PLON(:) = P__LON (1:ISIZE)
00625 ENDIF
00626 
00627 !  Allocate arrays to produce analysis increments  
00628 
00629 ALLOCATE (ZT2INC(ISIZE))
00630 ALLOCATE (ZH2INC(ISIZE))
00631 ALLOCATE (ZWSINC(ISIZE))
00632 ALLOCATE (ZWPINC(ISIZE))
00633 ALLOCATE (ZTLINC(ISIZE))
00634 ALLOCATE (ZTSINC(ISIZE))
00635 ALLOCATE (ZTPINC(ISIZE))
00636 ALLOCATE (ZSNINC(ISIZE))
00637 
00638 ! Screen-level innovations
00639 
00640 ZT2INC(:) = PT2M_O(:) - PTCLS(:)
00641 ZH2INC(:) = PHU2M_O(:) - PHCLS(:)
00642 
00643 ! Threshold for background check
00644 
00645 ZTHRES=RTHR_QC*SQRT(SIGWGO**2 + SIGWGB**2)
00646 
00647 ! Superficial soil moisture innovations in (m3/m3)
00648 
00649 DO J = 1, ISIZE
00650   IF (PWS_O(J) /= 999.0) THEN
00651     ZWGINC(J) = PWS_O(J) - PWS(J,1)
00652     IF (ABS(ZWGINC(J)) > ZTHRES) THEN 
00653       ZWGINC(J) = 0.0 ! background check
00654     ELSE
00655       INOBS = INOBS + 1
00656     ENDIF
00657   ELSE
00658     ZWGINC(J) = 0.0
00659   ENDIF
00660 ENDDO
00661 PRINT *,' NUMBER OF ASCAT OBSERVATIONS AFTER BACKGROUND CHECK  :: ',INOBS
00662 
00663 PRINT *,'           '
00664 PRINT *,'Mean T2m increments  ',SUM(ZT2INC)/SIZE(ZT2INC)
00665 PRINT *,'Mean HU2m increments ',SUM(ZH2INC)/SIZE(ZH2INC)
00666 PRINT *,'           '
00667 
00668 ! Interface (define arrays and perform unit conversions)
00669 
00670 PARG(:) = PARG(:)*100.0
00671 PSAB(:) = PSAB(:)*100.0
00672 
00673 ZWS(:) = XUNDEF
00674 ZWP(:) = XUNDEF
00675 ZTL(:) = XUNDEF
00676 
00677 WHERE (PWS(:,1)/=XUNDEF) 
00678   ZWS(:)      = PWS(:,1)*RD1*XRHOLW     ! conversion of m3/m3 -> mm
00679   ZWP(:)      = PWP(:,1)*PD2(:,1)*XRHOLW  ! conversion of m3/m3 -> mm
00680   ZTL(:)      = PTL(:,1)*PD2(:,1)*XRHOLW  ! conversion of m3/m3 -> mm
00681 END WHERE
00682 
00683 ZTCLS(:)    = PTCLS(:)
00684 ZHCLS(:)    = PHCLS(:)
00685 ZUCLS(:)    = PUCLS(:)
00686 ZVCLS(:)    = PVCLS(:)
00687 PSSTC(:)    = PTS_O(:)
00688 PWPINC1(:)  = XUNDEF
00689 PWPINC2(:)  = XUNDEF
00690 PWPINC3(:)  = XUNDEF
00691 PT2MBIAS(:) = XUNDEF
00692 PH2MBIAS(:) = XUNDEF
00693 
00694 ! Sea-ice surface properties
00695 
00696 PALBF(:)    = XUNDEF
00697 PEMISF(:)   = XUNDEF
00698 PZ0F(:)     = XUNDEF
00699 PZ0H(:)     = XUNDEF
00700 
00701 ! Climatological arrays set to missing values
00702 
00703 PSNC(:)     =  PSNS(:,1) ! need to read the snow climatology 
00704 PWSC(:)     =  XUNDEF
00705 PWPC(:)     =  XUNDEF
00706 PTSC(:)     =  XUNDEF
00707 PTPC(:)     =  XUNDEF
00708 
00709 DO J = 1, ISIZE
00710   PGELAT(J)   = PLAT(J) 
00711   PGELAM(J)   = PLON(J) 
00712   PGEMU(J)    = SIN(PLAT(J)*XPI/180.)
00713 ENDDO
00714 
00715 ZEVAP(:)   =  (PEVAP(:)/XLVTT*XDAY)/(NECHGU*3600.) ! conversion W/m2 -> mm/day
00716 ZEVAPTR(:) =  PEVAPTR(:)*XDAY
00717 ZSNS(:)    =  PSNS(:,1)
00718 
00719 DO J = 1, ISIZE
00720   ZTS(J) = PTS(J,1)
00721   ZTP(J) = PTP(J,1)
00722 ENDDO
00723 
00724 IDAT = IYEAR*10000. + IMONTH*100. + IDAY
00725 
00726 IF (LDINLINE) THEN
00727 ! Avoid division by zero in next WHERE statement; 
00728 ! this may occur in the extension zone
00729   WHERE (LD_MASKEXT (1:ISIZE))
00730     PD2(:,1) = 1.
00731   END WHERE
00732 ENDIF
00733 
00734 !  Soil analysis based on optimal interpolation
00735 
00736  CALL OI_CACSTS(ISIZE,ZT2INC,ZH2INC,ZWGINC,PWS_O,                       &
00737                 IDAT,NSSSSS,                                           &
00738                 ZTP,ZWP,ZTL,ZSNS,ZTS,ZWS,                              &
00739                 ZTCLS,ZHCLS,ZUCLS,ZVCLS,PSSTC,PWPINC1,PWPINC2,PWPINC3, &
00740                 PT2MBIAS,PH2MBIAS,                                     &
00741                 PRRCL,PRRSL,PRRCN,PRRSN,PATMNEB,ZEVAP,ZEVAPTR,         &
00742                 PITM,PVEG(:,1),PALBF,PEMISF,PZ0F,                      &
00743                 PIVEG,PARG,PD2(:,1),PSAB,PLAI(:,1),PRSMIN(:,1),PZ0H,   &
00744                 PTSC,PTPC,PWSC,PWPC,PSNC,                              &
00745                 PGELAT,PGELAM,PGEMU) 
00746 
00747 !  Store increments
00748 
00749 ZWSINC(:) = 0.0
00750 ZWPINC(:) = 0.0
00751 ZTLINC(:) = 0.0
00752 ZSNINC(:) = 0.0
00753 
00754 WHERE (PWS(:,1)/=XUNDEF)
00755   ZWSINC(:) = ZWS(:) - PWS(:,1)*(RD1*XRHOLW)    
00756   ZWPINC(:) = ZWP(:) - PWP(:,1)*(PD2(:,1)*XRHOLW) 
00757   ZTLINC(:) = ZTL(:) - PTL(:,1)*(PD2(:,1)*XRHOLW) 
00758   ZSNINC(:) = ZSNS(:) - PSNS(:,1)
00759 END WHERE
00760 
00761 IF (LDINLINE) THEN
00762 ! Avoid division by zero in next WHERE statement; 
00763 ! this may occur in the extension zone
00764   WHERE (LD_MASKEXT (1:ISIZE))
00765     PD2(:,1) = 1.
00766   END WHERE
00767 ENDIF
00768 
00769 !  Define soil moiture analyses over NATURE points
00770 
00771 WHERE (PWS(:,1)/=XUNDEF)
00772   PWS(:,1)  = ZWS(:)/(RD1*XRHOLW)
00773   PWP(:,1)  = ZWP(:)/(PD2(:,1)*XRHOLW)
00774   PTL(:,1)  = ZTL(:)/(PD2(:,1)*XRHOLW)
00775   PSNS(:,1) = ZSNS(:)
00776 END WHERE
00777 
00778 !  Perform temperature analysis according to surface types
00779 
00780 OINTERP_LST(:) = .FALSE.
00781 OINTERP_SST(:) = .FALSE.
00782 
00783 ZTSINC(:) = 0.0
00784 ZTPINC(:) = 0.0
00785 
00786 ! a) Temperature analysis of NATURE points
00787 
00788 WHERE (PTS(:,1)/=XUNDEF)
00789   ZTSINC(:) = ZTS(:) - PTS(:,1)
00790   ZTPINC(:) = ZTP(:) - PTP(:,1)
00791   PTS(:,1)  = ZTS(:)
00792   PTP(:,1)  = ZTP(:)
00793 END WHERE
00794 
00795 ! b) Temperature analysis of SEA and LAKE points
00796 
00797 DO J = 1, ISIZE
00798   IF (PITM(J) < 0.5) THEN
00799     IF (PSST(J)/=XUNDEF) THEN
00800       ZTSINC(J) = PTS_O(J) - PSST(J)
00801       PSST(J) = PTS_O(J)   ! canari
00802     ENDIF
00803     IF (PLST(J)/=XUNDEF) THEN
00804       PLST(J) = PTS_O(J)   ! canari
00805     ENDIF
00806   ELSE
00807     IF (PSST(J)/=XUNDEF) THEN
00808       PSST(J) = XUNDEF
00809       OINTERP_SST(J) = .TRUE.
00810     ENDIF
00811     IF (PLST(J)/=XUNDEF) THEN
00812       PLST(J) = XUNDEF
00813       OINTERP_LST(J) = .TRUE.
00814     ENDIF
00815   ENDIF
00816 ENDDO
00817 
00818 ! c) Temperature analysis of TOWN points
00819 
00820 WHERE (PTRD3(:)/=XUNDEF)
00821   PTRD3(:) = PTRD3(:) + ZT2INC(:)/(2.0*XPI)
00822 END WHERE
00823 
00824 ! Search for the nearest grid point values for lake and sea points
00825 ! at locations where the water fraction is less than 50 % 
00826 ! and therefore no useful information is given from the SST analysis
00827 ! A standard temperature gradient is applied to account for the atitude differences
00828 
00829 IF (LDINLINE) THEN
00830 
00831   IF (LLKEEPEXTZONE) THEN
00832 
00833     ZLST(:) = PLST(:)
00834      
00835     IF (LDINLINE) THEN
00836       WHERE (LD_MASKEXT (1:ISIZE))
00837         ZLST = XUNDEF
00838       END WHERE
00839     ENDIF
00840      
00841     CALL OI_HOR_EXTRAPOL_SURF(ISIZE,PLAT,PLON,ZLST,PLAT,PLON,PLST,OINTERP_LST,ZALT)
00842      
00843     ZSST(:) = PSST(:)
00844      
00845     IF (LDINLINE) THEN
00846       WHERE (LD_MASKEXT (1:ISIZE))
00847         ZSST = XUNDEF
00848       END WHERE
00849     ENDIF
00850      
00851     CALL OI_HOR_EXTRAPOL_SURF(ISIZE,PLAT,PLON,ZSST,PLAT,PLON,PSST,OINTERP_SST,ZALT)
00852 
00853   ELSE
00854 
00855     ISIZE1 = COUNT (.NOT. LD_MASKEXT)
00856    
00857     ALLOCATE (PSST1 (ISIZE1), PLST1 (ISIZE1), ZSST1 (ISIZE1), ZLST1 (ISIZE1), PLAT1 (ISIZE1), &
00858             & PLON1 (ISIZE1), ZALT1 (ISIZE1), OINTERP_LST1 (ISIZE1), OINTERP_SST1 (ISIZE1))
00859    
00860     ! remove extension zone
00861     J = 1
00862     DO J1 = 1, ISIZE
00863       IF (.NOT. LD_MASKEXT (J1)) THEN
00864         PSST1 (J) = PSST (J1)
00865         PLST1 (J) = PLST (J1)
00866         PLAT1 (J) = PLAT (J1)
00867         PLON1 (J) = PLON (J1)
00868         ZALT1 (J) = ZALT (J1)
00869         OINTERP_LST1 (J) = OINTERP_LST (J1)
00870         OINTERP_SST1 (J) = OINTERP_SST (J1)
00871         J = J + 1
00872       ENDIF
00873     ENDDO
00874    
00875     ZLST1(:) = PLST1(:)
00876     CALL OI_HOR_EXTRAPOL_SURF(ISIZE1,PLAT1,PLON1,ZLST1,PLAT1,PLON1,PLST1,OINTERP_LST1,ZALT1)
00877      
00878     ZSST1(:) = PSST1(:)
00879     CALL OI_HOR_EXTRAPOL_SURF(ISIZE1,PLAT1,PLON1,ZSST1,PLAT1,PLON1,PSST1,OINTERP_SST1,ZALT1)
00880    
00881     ! copy back
00882     J = 1
00883     DO J1 = 1, ISIZE
00884       IF (.NOT. LD_MASKEXT (J1)) THEN
00885         PSST (J1) = PSST1 (J)
00886         PLST (J1) = PLST1 (J) 
00887         J = J + 1
00888       ENDIF
00889     ENDDO
00890    
00891     DEALLOCATE (PSST1, PLST1, ZSST1, ZLST1, PLAT1, PLON1, &
00892               & ZALT1, OINTERP_LST1, OINTERP_SST1)
00893 
00894   ENDIF
00895 
00896 ELSE
00897  
00898   ZLST(:) = PLST(:)
00899   CALL OI_HOR_EXTRAPOL_SURF(ISIZE,PLAT,PLON,ZLST,PLAT,PLON,PLST,OINTERP_LST,ZALT)
00900 
00901   ZSST(:) = PSST(:)
00902   CALL OI_HOR_EXTRAPOL_SURF(ISIZE,PLAT,PLON,ZSST,PLAT,PLON,PSST,OINTERP_SST,ZALT)
00903 
00904 ENDIF
00905 
00906 ! PRINT values produced by OI_HO_EXTRAPOL_SURF
00907 
00908 IF (LPRINT) THEN
00909   DO J = 1, ISIZE
00910     IF (OINTERP_LST(J)) THEN
00911       PRINT *,'Lake surface temperature set to ',PLST(J),'from nearest neighbour at J=',J
00912     ENDIF
00913     IF (OINTERP_SST(J)) THEN
00914       PRINT *,'Sea surface temperature set to ',PSST(J),'from nearest neighbour at J=',J
00915     ENDIF
00916   ENDDO
00917 ENDIF
00918 
00919 ! PRINT statistics of the soil analysis
00920 
00921 PRINT *,'---------------------------------------------------------------'
00922 PRINT *,'Mean WS increments over NATURE ',SUM(ZWSINC,XNATURE > 0.)/NDIM_NATURE
00923 PRINT *,'Mean WP increments over NATURE ',SUM(ZWPINC,XNATURE > 0.)/NDIM_NATURE
00924 PRINT *,'Mean TS increments over NATURE ',SUM(ZTSINC,XNATURE > 0.)/NDIM_NATURE
00925 PRINT *,'Mean TP increments over NATURE ',SUM(ZTPINC,XNATURE > 0.)/NDIM_NATURE
00926 PRINT *,'Mean TL increments over NATURE ',SUM(ZTLINC,XNATURE > 0.)/NDIM_NATURE
00927 PRINT *,'Mean SN increments over NATURE ',SUM(ZSNINC,XNATURE > 0.)/NDIM_NATURE
00928 PRINT *,'---------------------------------------------------------------'
00929 PRINT *,'Mean WS increments over SEA    ',SUM(ZWSINC,XSEA > 0.)/NDIM_SEA
00930 PRINT *,'Mean WP increments over SEA    ',SUM(ZWPINC,XSEA > 0.)/NDIM_SEA
00931 PRINT *,'Mean TS increments over SEA    ',SUM(ZTSINC,XSEA > 0.)/NDIM_SEA
00932 PRINT *,'Mean TP increments over SEA    ',SUM(ZTPINC,XSEA > 0.)/NDIM_SEA
00933 PRINT *,'Mean TL increments over SEA    ',SUM(ZTLINC,XSEA > 0.)/NDIM_SEA
00934 PRINT *,'Mean SN increments over SEA    ',SUM(ZSNINC,XSEA > 0.)/NDIM_SEA
00935 PRINT *,'---------------------------------------------------------------'
00936 
00937 IF (.NOT. LDINLINE) THEN 
00938 !   Write analysis in LFI file PREP
00939 #ifdef LFI
00940   CFILEOUT_LFI='PREP'
00941 #endif
00942 ENDIF
00943 
00944  CALL FLAG_UPDATE(.FALSE.,.FALSE.,.TRUE.,.FALSE.)
00945  CALL INIT_IO_SURF_n(YPROGRAM,'FULL  ','SURF  ','WRITE')
00946 
00947 IF (LDINLINE) THEN
00948 #ifdef ARO
00949 ! Count 2D fields in MSE
00950   NCOUNTW_TOT = 0
00951   LWRITE      = .FALSE.
00952   LCOUNTW     = .TRUE.
00953   NCOUNTW     = 0
00954 ! NINDX1, NINDX2, NKPROMA already set 
00955   IF (.NOT. LFMWRIT) &
00956  & CALL SURFEX_FIELD_BUF_SET_RECORD (YSURFEX_CACHE_OUT, .FALSE.)
00957 
00958   CALL WRITE
00959    
00960   IF (LFMWRIT) THEN
00961     ALLOCATE (XGPGW (IGPCOMP, NCOUNTW_TOT))
00962     XGPGW = XUNDEF
00963   ELSE
00964     CALL SURFEX_FIELD_BUF_PREALLOC (YSURFEX_CACHE_OUT)
00965     CALL SURFEX_FIELD_BUF_SET_RECORD (YSURFEX_CACHE_OUT, .TRUE.)
00966   ENDIF
00967 
00968   LWRITE      = .TRUE.
00969   LCOUNTW     = .FALSE.
00970   NCOUNTW     = 0
00971 #endif
00972 
00973 ENDIF
00974 
00975  CALL WRITE
00976 
00977 IF (LDINLINE) THEN
00978 #ifdef ARO
00979   IF (LFMWRIT) DEALLOCATE (XGPGW)
00980 #endif
00981 ENDIF
00982 
00983  CALL END_IO_SURF_n(YPROGRAM)
00984  CALL IO_BUFF_CLEAN_n
00985 
00986 DEALLOCATE(NWORK)
00987 DEALLOCATE(XWORK)
00988 DEALLOCATE(XWORK2)
00989 
00990 IF (.NOT. LDINLINE) THEN
00991   PRINT *,'after write in PREP file'
00992 ENDIF
00993 
00994 ! -------------------------------------------------------------------------------------
00995 IF (LHOOK) CALL DR_HOOK ('OI_CONTROL', 1, ZHOOK_HANDLE)
00996 
00997 CONTAINS
00998 
00999 SUBROUTINE WRITE 
01000 
01001 REAL(KIND=JPRB) :: ZHOOK_HANDLE
01002 
01003 IF (LHOOK) CALL DR_HOOK ('OI_CONTROL:WRITE', 0, ZHOOK_HANDLE)
01004 
01005  CALL DD ('WG1', PWS (:,1))
01006 
01007 YVAR='WG1'
01008 YPREFIX='X_Y_WG1 (m3/m3)                                   '
01009  CALL WRITE_SURF(YPROGRAM,YVAR,PWS,IRESP,HCOMMENT=YPREFIX)
01010 
01011  CALL DD ('WG2', PWP (:,1))
01012 
01013 YVAR='WG2'
01014 YPREFIX='X_Y_WG2 (m3/m3)                                   '
01015  CALL WRITE_SURF(YPROGRAM,YVAR,PWP,IRESP,HCOMMENT=YPREFIX)
01016 
01017  CALL DD ('WGI2', PTL (:,1))
01018 
01019 YVAR='WGI2'
01020 YPREFIX='X_Y_WGI2 (m3/m3)                                  '
01021  CALL WRITE_SURF(YPROGRAM,YVAR,PTL,IRESP,HCOMMENT=YPREFIX)
01022 
01023  CALL DD ('TG1', PTS (:,1))
01024 
01025 YVAR='TG1'
01026 YPREFIX='X_Y_TG1 (K)                                       '
01027  CALL WRITE_SURF(YPROGRAM,YVAR,PTS,IRESP,HCOMMENT=YPREFIX)
01028 
01029  CALL DD ('TG2', PTP (:,1))
01030 
01031 YVAR='TG2'
01032 YPREFIX='X_Y_TG2 (K)                                       '
01033  CALL WRITE_SURF(YPROGRAM,YVAR,PTP,IRESP,HCOMMENT=YPREFIX)
01034 
01035  CALL DD ('SST', PSST)
01036 
01037 YVAR='SST'
01038 YPREFIX='X_Y_SST (K)                                       '
01039  CALL WRITE_SURF(YPROGRAM,YVAR,PSST,IRESP,HCOMMENT=YPREFIX)
01040 
01041  CALL DD ('TS_WATER', PLST)
01042 
01043 YVAR='TS_WATER'
01044 YPREFIX='X_Y_TS_WATER (K)                                  '
01045  CALL WRITE_SURF(YPROGRAM,YVAR,PLST,IRESP,HCOMMENT=YPREFIX)
01046 
01047 IF (NSIZE_TOWN > 0 .AND. LAROME) THEN
01048   CALL DD ('T_ROAD3', PTRD3)
01049 
01050   YVAR='TROAD3'
01051   YPREFIX='X_Y_T_ROAD3 (K)                                   '
01052   CALL WRITE_SURF(YPROGRAM,YVAR,PTRD3,IRESP,HCOMMENT=YPREFIX)
01053 ENDIF
01054 
01055  CALL DD ('WSNOW_VEG1', PSNS (:,1))
01056 
01057 YVAR='WSN_VEG1'
01058 YPREFIX='X_Y_WSNOW_VEG1 (kg/m2)                            '
01059  CALL WRITE_SURF(YPROGRAM,YVAR,PSNS,IRESP,HCOMMENT=YPREFIX)
01060 
01061 IF (LHOOK) CALL DR_HOOK ('OI_CONTROL:WRITE', 1, ZHOOK_HANDLE)
01062 
01063 END SUBROUTINE WRITE
01064 
01065 SUBROUTINE DD (CDN, PX)
01066  CHARACTER(LEN=*), INTENT (IN) :: CDN
01067 REAL, INTENT (IN) :: PX (:)
01068 
01069 REAL :: ZX (SIZE (PX))
01070 INTEGER :: I, N
01071 REAL(KIND=JPRB) :: ZHOOK_HANDLE
01072 
01073 IF (LHOOK) CALL DR_HOOK ('OI_CONTROL:DD', 0, ZHOOK_HANDLE)
01074 
01075 IF (LDINLINE) THEN
01076 #ifdef ARO
01077   IF (.NOT.LWRITE.AND.LHOOK) CALL DR_HOOK ('OI_CONTROL:DD', 1, ZHOOK_HANDLE)
01078   IF (.NOT.LWRITE) RETURN
01079 #endif
01080   N = COUNT (.NOT. LD_MASKEXT)
01081   ZX (1:N) = PACK (PX, .NOT. LD_MASKEXT)
01082 ELSE
01083   ZX = PX
01084   N = SIZE (PX)
01085 ENDIF
01086 
01087 WRITE (0, *) TRIM(CDN)//" = " 
01088 WRITE (0, *) N, MINVAL(ZX(1:N)), MAXVAL(ZX(1:N))
01089 !WRITE (0, '(10(E14.6,", "))') ZX (1:N)
01090 
01091 IF (LHOOK) CALL DR_HOOK ('OI_CONTROL:DD', 1, ZHOOK_HANDLE)
01092 
01093 END SUBROUTINE DD
01094 
01095 END SUBROUTINE OI_CONTROL