Juan 24/05/2017: add MPPDB_CHECK on SURFEX routine for reprod check with MNH_PARALLEL key
[MNH-git_open_source-lfs.git] / src / SURFEX / zoom_pgd_cover.F90
1 !SURFEX_LIC Copyright 1994-2014 Meteo-France 
2 !SURFEX_LIC This is part of the SURFEX software governed by the CeCILL-C  licence
3 !SURFEX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
4 !SURFEX_LIC for details. version 1.
5 !     #########
6       SUBROUTINE ZOOM_PGD_COVER(HPROGRAM,HINIFILE,HINIFILETYPE,OECOCLIMAP)
7 !     ###########################################################
8
9 !!
10 !!    PURPOSE
11 !!    -------
12 !!   This program prepares the physiographic data fields.
13 !!
14 !!    METHOD
15 !!    ------
16 !!   
17 !!    EXTERNAL
18 !!    --------
19 !!
20 !!
21 !!    IMPLICIT ARGUMENTS
22 !!    ------------------
23 !!
24 !!
25 !!    REFERENCE
26 !!    ---------
27 !!
28 !!    AUTHOR
29 !!    ------
30 !!
31 !!    V. Masson                   Meteo-France
32 !!
33 !!    MODIFICATION
34 !!    ------------
35 !!
36 !!    Original     13/10/03
37 !     Modification 17/04/12 M.Tomasini All COVER physiographic fields are now 
38 !!                                     interpolated for spawning => 
39 !!                                     ABOR1_SFX if (.NOT.OECOCLIMAP) in comment
40 !     Modification 05/02/15 M.Moge : MPPDB_CHECK + use NSIZE_FULL instead of SIZE(XLAT) (for clarity)
41 !!      J.Escobar 18/12/2015 : missing interface
42 !----------------------------------------------------------------------------
43 !
44 !*    0.     DECLARATION
45 !            -----------
46 !
47 USE MODD_DATA_COVER_PAR,   ONLY : JPCOVER
48 USE MODD_SURF_ATM_GRID_n,  ONLY : XLAT, XLON, CGRID, XGRID_PAR
49 USE MODD_SURF_ATM_n,       ONLY : XCOVER, LCOVER, XSEA, XWATER, XNATURE, XTOWN, &
50                                     NSIZE_NATURE, NSIZE_SEA, NR_NATURE, NR_SEA, &
51                                     NSIZE_TOWN, NSIZE_WATER,NR_TOWN,NR_WATER,NSIZE_FULL,&
52                                     NDIM_NATURE, NDIM_SEA,                  &
53                                     NDIM_TOWN,NDIM_WATER,NDIM_FULL  
54 USE MODD_PREP,             ONLY : CINGRID_TYPE, CINTERP_TYPE
55 !
56 USE MODI_CONVERT_COVER_FRAC
57 USE MODI_OPEN_AUX_IO_SURF
58 USE MODI_READ_SURF
59 USE MODI_CLOSE_AUX_IO_SURF
60 USE MODI_PREP_GRID_EXTERN
61 USE MODI_HOR_INTERPOL
62 USE MODI_HOR_INTERPOL_1COV
63 USE MODI_PREP_OUTPUT_GRID
64 USE MODI_OLD_NAME
65 USE MODI_SUM_ON_ALL_PROCS
66 USE MODI_GET_LUOUT
67 USE MODI_CLEAN_PREP_OUTPUT_GRID
68 USE MODI_GET_1D_MASK
69 USE MODI_READ_LCOVER
70 USE MODI_READ_SURFX2COV_1COV_MNH
71 !
72 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
73 USE PARKIND1  ,ONLY : JPRB
74 !
75 #ifdef MNH_PARALLEL
76 USE MODE_MPPDB
77 !
78 #endif
79 IMPLICIT NONE
80 !
81 !*    0.1    Declaration of dummy arguments
82 !            ------------------------------
83 !
84  CHARACTER(LEN=6),     INTENT(IN)  :: HPROGRAM    ! program calling
85  CHARACTER(LEN=28),    INTENT(IN)  :: HINIFILE    ! input atmospheric file name
86  CHARACTER(LEN=6),     INTENT(IN)  :: HINIFILETYPE! input atmospheric file type
87 LOGICAL,              INTENT(OUT) :: OECOCLIMAP  ! flag to use ecoclimap
88 !
89 !
90 !*    0.2    Declaration of local variables
91 !            ------------------------------
92 !
93 INTEGER :: IRESP
94 INTEGER :: ILUOUT
95 INTEGER :: INI     ! total 1D dimension (input grid)
96 INTEGER :: IL      ! total 1D dimension (output grid)
97 INTEGER :: JCOVER  ! loop counter
98 INTEGER :: IVERSION       ! surface version
99 #ifdef MNH_PARALLEL
100 REAL, DIMENSION(:), POINTER     :: ZCOVER
101 #else
102 REAL, DIMENSION(:,:), POINTER     :: ZCOVER
103 #endif
104 REAL, DIMENSION(:,:), POINTER :: ZSEA1, ZWATER1, ZNATURE1, ZTOWN1
105 REAL, DIMENSION(:,:), POINTER :: ZSEA2, ZWATER2, ZNATURE2, ZTOWN2
106 REAL, DIMENSION(:),   ALLOCATABLE :: ZSUM
107 CHARACTER(LEN=16) :: YRECFM         ! Name of the article to be read
108 CHARACTER(LEN=100) :: YCOMMENT
109 REAL(KIND=JPRB) :: ZHOOK_HANDLE
110 !------------------------------------------------------------------------------
111 IF (LHOOK) CALL DR_HOOK('ZOOM_PGD_COVER',0,ZHOOK_HANDLE)
112  CALL GET_LUOUT(HPROGRAM,ILUOUT)
113 !
114 !*      1.     Preparation of IO for reading in the file
115 !              -----------------------------------------
116 !
117 !* Note that all points are read, even those without physical meaning.
118 !  These points will not be used during the horizontal interpolation step.
119 !  Their value must be defined as XUNDEF.
120 !
121  CALL OPEN_AUX_IO_SURF(HINIFILE,HINIFILETYPE,'FULL  ')
122 !
123  CALL READ_SURF(HPROGRAM,'ECOCLIMAP',OECOCLIMAP,IRESP)
124 !
125 !------------------------------------------------------------------------------
126 !
127 !*      2.     Reading of grid
128 !              ---------------
129 !
130  CALL PREP_GRID_EXTERN(HINIFILETYPE,ILUOUT,CINGRID_TYPE,CINTERP_TYPE,INI)
131 !
132  CALL PREP_OUTPUT_GRID(ILUOUT,CGRID,XGRID_PAR,XLAT,XLON)
133 #ifdef MNH_PARALLEL
134  CALL MPPDB_CHECK_SURFEX2D(XLAT,"ZOOM_PGD_COVER:XLAT",PRECISION,ILUOUT)
135  CALL MPPDB_CHECK_SURFEX2D(XLON,"ZOOM_PGD_COVER:XLON",PRECISION,ILUOUT)
136 #endif
137 !
138 !------------------------------------------------------------------------------
139 !
140 !*      3.     Reading of cover
141 !              ----------------
142 !
143 YRECFM='VERSION'
144 CALL READ_SURF(HPROGRAM,YRECFM,IVERSION,IRESP)
145 !
146 ALLOCATE(LCOVER(JPCOVER))
147 ALLOCATE(ZSEA1   (INI,1))
148 ALLOCATE(ZNATURE1(INI,1))
149 ALLOCATE(ZWATER1 (INI,1))
150 ALLOCATE(ZTOWN1  (INI,1))
151 !
152 IF (IVERSION>=7) THEN
153   CALL READ_SURF(HPROGRAM,'FRAC_SEA   ',ZSEA1(:,1),   IRESP,HDIR='A')
154   CALL READ_SURF(HPROGRAM,'FRAC_NATURE',ZNATURE1(:,1),IRESP,HDIR='A')
155   CALL READ_SURF(HPROGRAM,'FRAC_WATER ',ZWATER1(:,1), IRESP,HDIR='A')
156   CALL READ_SURF(HPROGRAM,'FRAC_TOWN  ',ZTOWN1(:,1),  IRESP,HDIR='A')
157   !
158   CALL OLD_NAME(HPROGRAM,'COVER_LIST      ',YRECFM)
159 !  CALL READ_SURF(HPROGRAM,YRECFM,LCOVER(:),IRESP,HDIR='-')
160   CALL READ_LCOVER(HPROGRAM,LCOVER)
161   !
162 #ifdef MNH_PARALLEL
163   ALLOCATE(ZCOVER(INI))
164 #else
165   ALLOCATE(ZCOVER(INI,JPCOVER))
166 #endif
167   !
168 ELSE
169 #ifdef MNH_PARALLEL
170   ! we assume that IVERSION>=7
171 #else
172   CALL OLD_NAME(HPROGRAM,'COVER_LIST      ',YRECFM)
173 !  CALL READ_SURF(HPROGRAM,YRECFM,LCOVER(:),IRESP,HDIR='-')
174   CALL READ_LCOVER(HPROGRAM,LCOVER)
175   !
176   ALLOCATE(ZCOVER(INI,JPCOVER))
177   CALL READ_SURF(HPROGRAM,YRECFM,ZCOVER(:,:),LCOVER,IRESP,HDIR='A')
178   !
179   CALL CONVERT_COVER_FRAC(ZCOVER,ZSEA1(:,1),ZNATURE1(:,1),ZTOWN1(:,1),ZWATER1(:,1))
180 #endif
181 ENDIF
182 !
183 !------------------------------------------------------------------------------
184 !
185 !*      4.     Reading of cover & Interpolations
186 !              --------------
187 !
188 IL = NSIZE_FULL
189 ALLOCATE(XCOVER(IL,JPCOVER))
190 ALLOCATE(ZSUM(IL))
191 ZSUM = 0.
192 !
193 ! on lit les cover une apres l'autre, et on appelle hor_interpol sur chaque cover separement
194 !
195 #ifdef MNH_PARALLEL
196 IF ( HPROGRAM == 'MESONH' ) THEN
197   DO JCOVER=1,JPCOVER
198     IF ( LCOVER( JCOVER ) ) THEN
199       CALL READ_SURFX2COV_1COV_MNH(YRECFM,INI,JCOVER,ZCOVER(:),IRESP,YCOMMENT,'A')
200     ELSE
201       ZCOVER(:) = 0.
202     ENDIF
203     !
204     CALL HOR_INTERPOL_1COV(ILUOUT,ZCOVER,XCOVER(:,JCOVER))
205 #ifdef MNH_PARALLEL
206     CALL MPPDB_CHECK_SURFEX3D(XCOVER,"ZOOM_PGD_COVER:XCOVER",PRECISION,ILUOUT,'FULL',JPCOVER)
207 #endif
208     !
209     ZSUM(:) = ZSUM(:) + XCOVER(:,JCOVER)
210     !
211   ENDDO
212 ELSE
213   
214 ENDIF
215 #else
216  CALL HOR_INTERPOL(ILUOUT,ZCOVER,XCOVER)
217 #endif
218 !
219 !  Coherence check
220 !
221 #ifdef MNH_PARALLEL
222 CALL MPPDB_CHECK_SURFEX2D(ZSUM,"ZOOM_PGD_COVER:ZSUM",PRECISION,ILUOUT)
223 #endif
224 DO JCOVER=1,JPCOVER
225   XCOVER(:,JCOVER) = XCOVER(:,JCOVER)/ZSUM(:)
226   IF (ALL(XCOVER(:,JCOVER)==0.)) LCOVER(JCOVER) = .FALSE.
227 END DO
228 !
229 CALL CLOSE_AUX_IO_SURF(HINIFILE,HINIFILETYPE)
230 !
231 !
232 DEALLOCATE(ZCOVER)
233 !
234 ALLOCATE(ZSEA2  (IL,1))
235 ALLOCATE(ZNATURE2(IL,1))
236 ALLOCATE(ZWATER2 (IL,1))
237 ALLOCATE(ZTOWN2  (IL,1))
238 !
239  CALL HOR_INTERPOL(ILUOUT,ZSEA1,ZSEA2)
240  CALL HOR_INTERPOL(ILUOUT,ZNATURE1,ZNATURE2)
241  CALL HOR_INTERPOL(ILUOUT,ZWATER1,ZWATER2)
242  CALL HOR_INTERPOL(ILUOUT,ZTOWN1,ZTOWN2)
243 !
244 DEALLOCATE(ZSEA1)
245 DEALLOCATE(ZNATURE1)
246 DEALLOCATE(ZWATER1)
247 DEALLOCATE(ZTOWN1)
248 !
249 ALLOCATE(XSEA   (IL))
250 ALLOCATE(XNATURE(IL))
251 ALLOCATE(XWATER (IL))
252 ALLOCATE(XTOWN  (IL))
253 !
254 XSEA(:)   = ZSEA2   (:,1)
255 XNATURE(:)= ZNATURE2(:,1)
256 XWATER(:) = ZWATER2 (:,1)
257 XTOWN(:)  = ZTOWN2  (:,1)
258 !
259 DEALLOCATE(ZSEA2)
260 DEALLOCATE(ZNATURE2)
261 DEALLOCATE(ZWATER2)
262 DEALLOCATE(ZTOWN2)
263 !
264  CALL CLEAN_PREP_OUTPUT_GRID
265 !------------------------------------------------------------------------------
266 !
267 !*      6.     Fractions
268 !              ---------
269 !
270 ! When the model runs in multiproc, NSIZE* represents the number of points
271 ! on a proc, and NDIM* the total number of points on all procs.
272 ! The following definition of NDIM* won't be correct any more when the PGD
273 ! runs in multiproc.
274 !
275 NSIZE_NATURE    = COUNT(XNATURE(:) > 0.0)
276 NSIZE_WATER     = COUNT(XWATER (:) > 0.0)
277 NSIZE_SEA       = COUNT(XSEA   (:) > 0.0)
278 NSIZE_TOWN      = COUNT(XTOWN  (:) > 0.0)
279 NSIZE_FULL      = IL
280 !
281 NDIM_NATURE    = SUM_ON_ALL_PROCS(HPROGRAM,CGRID,XNATURE(:) > 0., 'DIM')
282 NDIM_WATER     = SUM_ON_ALL_PROCS(HPROGRAM,CGRID,XWATER (:) > 0., 'DIM')
283 NDIM_SEA       = SUM_ON_ALL_PROCS(HPROGRAM,CGRID,XSEA   (:) > 0., 'DIM')
284 NDIM_TOWN      = SUM_ON_ALL_PROCS(HPROGRAM,CGRID,XTOWN  (:) > 0., 'DIM')
285 ZSUM=1.
286 NDIM_FULL      = SUM_ON_ALL_PROCS(HPROGRAM,CGRID,ZSUM   (:) ==1., 'DIM')
287 DEALLOCATE(ZSUM)
288 !
289 ALLOCATE(NR_NATURE (NSIZE_NATURE))
290 ALLOCATE(NR_TOWN   (NSIZE_TOWN  ))
291 ALLOCATE(NR_WATER  (NSIZE_WATER ))
292 ALLOCATE(NR_SEA    (NSIZE_SEA   ))
293 !
294 IF (NSIZE_SEA   >0)CALL GET_1D_MASK( NSIZE_SEA,    NSIZE_FULL, XSEA   , NR_SEA   )
295 IF (NSIZE_WATER >0)CALL GET_1D_MASK( NSIZE_WATER,  NSIZE_FULL, XWATER , NR_WATER )
296 IF (NSIZE_TOWN  >0)CALL GET_1D_MASK( NSIZE_TOWN,   NSIZE_FULL, XTOWN  , NR_TOWN  )
297 IF (NSIZE_NATURE>0)CALL GET_1D_MASK( NSIZE_NATURE, NSIZE_FULL, XNATURE, NR_NATURE)
298 #ifdef MNH_PARALLEL
299 CALL MPPDB_CHECK_SURFEX2D(XSEA,"ZOOM_PGD_COVER:XSEA",PRECISION,ILUOUT)
300 CALL MPPDB_CHECK_SURFEX2D(XWATER,"ZOOM_PGD_COVER:XWATER",PRECISION,ILUOUT)
301 CALL MPPDB_CHECK_SURFEX2D(XTOWN,"ZOOM_PGD_COVER:XTOWN",PRECISION,ILUOUT)
302 CALL MPPDB_CHECK_SURFEX2D(XNATURE,"ZOOM_PGD_COVER:XNATURE",PRECISION,ILUOUT)
303 #endif
304 IF (LHOOK) CALL DR_HOOK('ZOOM_PGD_COVER',1,ZHOOK_HANDLE)
305
306 !_______________________________________________________________________________
307 !
308 END SUBROUTINE ZOOM_PGD_COVER