Juan 8/12/2016: add management of LEN_HREC in MNH & SURFEX
[MNH-git_open_source-lfs.git] / src / SURFEX / write_diag_pgd_isban.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 WRITE_DIAG_PGD_ISBA_n(HPROGRAM)
7 !     #########################################
8 !
9 !!****  *WRITE_DIAG_PGD_ISBA_n* - writes the ISBA physiographic diagnostic fields
10 !!
11 !!    PURPOSE
12 !!    -------
13 !!
14 !!**  METHOD
15 !!    ------
16 !!
17 !!    EXTERNAL
18 !!    --------
19 !!
20 !!
21 !!    IMPLICIT ARGUMENTS
22 !!    ------------------
23 !!
24 !!    REFERENCE
25 !!    ---------
26 !!
27 !!
28 !!    AUTHOR
29 !!    ------
30 !!      V. Masson   *Meteo France*      
31 !!
32 !!    MODIFICATIONS
33 !!    -------------
34 !!      Original    01/2004 
35 !!      Modified    10/2004 by P. Le Moigne: add XZ0REL, XVEGTYPE_PATCH
36 !!      Modified    11/2005 by P. Le Moigne: limit length of VEGTYPE_PATCH field names
37 !!      M.Moge    01/2016  using WRITE_SURF_FIELD2D/3D for 2D/3D surfex fields writes
38 !-------------------------------------------------------------------------------
39 !
40 !*       0.    DECLARATIONS
41 !              ------------
42 !
43 USE MODD_SURF_PAR,   ONLY : XUNDEF, NUNDEF
44 USE MODD_ISBA_n,     ONLY : NPATCH, CPHOTO, CHORT, CISBA,                           &
45                               XLAI, XVEG, XZ0,XALBNIR_SOIL,XALBVIS_SOIL,XALBUV_SOIL,&
46                               XRSMIN, XGAMMA, XRGL, XCV, XEMIS, XDG, XWRMAX_CF,     &
47                               XZ0REL, XVEGTYPE_PATCH, XALBNIR, XALBVIS, XALBUV,     &
48                               XPATCH, XWATSUP, TSEED, TREAP, XIRRIG, XD_ICE,        &
49                               XROOTFRAC, NWG_LAYER, XDROOT, XDG2,                   &
50                               XWSAT, XWFC, XWWILT, XRUNOFFD, CSOC, XFRACSOC   
51 USE MODD_AGRI,       ONLY : LAGRIP
52 !
53 USE MODD_DIAG_MISC_ISBA_n,ONLY : LSURF_DIAG_ALBEDO
54 !
55 USE MODD_IO_SURF_FA, ONLY : LFANOCOMPACT, LPREP
56 !
57 USE MODD_CH_ISBA_n,  ONLY : XSOILRC_SO2, XSOILRC_O3, CCH_DRY_DEP, NBEQ
58 USE MODI_INIT_IO_SURF_n
59 USE MODI_WRITE_SURF
60 USE MODI_WRITE_SURF_FIELD2D
61 USE MODI_END_IO_SURF_n
62 !
63 !
64 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
65 USE PARKIND1  ,ONLY : JPRB
66 !
67 IMPLICIT NONE
68 !
69 !*       0.1   Declarations of arguments
70 !              -------------------------
71 !
72  CHARACTER(LEN=6),  INTENT(IN)  :: HPROGRAM ! program calling
73 !
74 !*       0.2   Declarations of local variables
75 !              -------------------------------
76 !
77 REAL, DIMENSION(SIZE(XDG,1),SIZE(XDG,3)) :: ZWORK ! Work array
78 !
79 INTEGER           :: IRESP          ! IRESP  : return-code if a problem appears
80  CHARACTER(LEN=LEN_HREC) :: YRECFM         ! Name of the article to be read
81  CHARACTER(LEN=100):: YCOMMENT       ! Comment string
82  CHARACTER(LEN=100):: YCOMMENTUNIT   ! Comment string : unit of the datas in the field to write
83  CHARACTER(LEN=2)  :: YLVLV, YPAS
84 !
85 INTEGER           :: JJ, JL, JP, ILAYER
86 REAL(KIND=JPRB) :: ZHOOK_HANDLE
87 !-------------------------------------------------------------------------------
88 !
89 !         Initialisation for IO
90 !
91 IF (LHOOK) CALL DR_HOOK('WRITE_DIAG_PGD_ISBA_N',0,ZHOOK_HANDLE)
92  CALL INIT_IO_SURF_n(HPROGRAM,'NATURE','ISBA  ','WRITE')
93 !
94 !-------------------------------------------------------------------------------
95 !
96 !* Leaf Area Index
97 !
98 IF (CPHOTO=='NON' .OR. CPHOTO=='AGS' .OR. CPHOTO=='AST') THEN
99   !
100   YRECFM='LAI'
101   YCOMMENT='leaf area index'
102   YCOMMENTUNIT='-'
103   !
104   CALL WRITE_SURF_FIELD2D(HPROGRAM,XLAI(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
105   !
106 ENDIF
107 !
108 !-------------------------------------------------------------------------------
109 !
110 !* Vegetation fraction
111 !
112 YRECFM='VEG'
113 YCOMMENT='vegetation fraction'
114 YCOMMENTUNIT='-'
115 !
116 CALL WRITE_SURF_FIELD2D(HPROGRAM,XVEG(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
117 !
118 !* Surface roughness length (without snow)
119 !
120 YRECFM='Z0VEG'
121 YCOMMENT='surface roughness length (without snow)'
122 YCOMMENTUNIT='M'
123 !
124 CALL WRITE_SURF_FIELD2D(HPROGRAM,XZ0(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
125 !
126 !-------------------------------------------------------------------------------
127 !
128 !* Fraction for each patch
129 !
130 IF(.NOT.LFANOCOMPACT.OR.LPREP)THEN
131   YRECFM='PATCH'
132   YCOMMENT='fraction for each patch'
133   YCOMMENTUNIT='-'
134   CALL WRITE_SURF_FIELD2D(HPROGRAM,XPATCH(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
135 ENDIF
136 !-------------------------------------------------------------------------------
137 !
138 !* Soil depth for each patch
139 !
140 DO JL=1,SIZE(XDG,2)
141   IF (JL<10) THEN
142     WRITE(YRECFM,FMT='(A2,I1)') 'DG',JL
143   ELSE
144     WRITE(YRECFM,FMT='(A2,I2)') 'DG',JL          
145   ENDIF
146   YCOMMENT='soil depth'
147   YCOMMENTUNIT='M'
148   CALL WRITE_SURF_FIELD2D(HPROGRAM,XDG(:,JL,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
149 END DO
150 !-------------------------------------------------------------------------------
151 !
152 IF(CISBA=='DIF')THEN
153 !
154 !* Root depth
155 !
156   YRECFM='DROOT_DIF'
157   YCOMMENT='Root depth in ISBA-DIF'
158 !
159   YCOMMENTUNIT='-'
160   CALL WRITE_SURF_FIELD2D(HPROGRAM,XDROOT(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
161 !
162   YRECFM='DG2_DIF'
163   YCOMMENT='DG2 depth in ISBA-DIF'
164 !
165   YCOMMENTUNIT='-'
166   CALL WRITE_SURF_FIELD2D(HPROGRAM,XDG2(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
167 !
168 !* Runoff depth
169 !
170   YRECFM='RUNOFFD'
171   YCOMMENT='Runoff deph in ISBA-DIF'
172 !
173   YCOMMENTUNIT='-'
174   CALL WRITE_SURF_FIELD2D(HPROGRAM,XRUNOFFD(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
175 !
176 !* Total soil depth for mositure
177 !
178   ZWORK(:,:)=XUNDEF
179   DO JP=1,SIZE(XDG,3)
180      DO JJ=1,SIZE(XDG,1)
181         JL=NWG_LAYER(JJ,JP)
182         IF(JL/=NUNDEF)THEN
183           ZWORK(JJ,JP)=XDG(JJ,JL,JP)
184         ENDIF
185      ENDDO
186   ENDDO
187   YRECFM='DTOT_DIF'
188   YCOMMENT='Total soil depth for moisture in ISBA-DIF'
189   YCOMMENTUNIT='-'
190   CALL WRITE_SURF_FIELD2D(HPROGRAM,ZWORK(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
191 !
192 !* Root fraction for each patch
193 !
194   DO JL=1,SIZE(XROOTFRAC,2)
195      IF (JL<10) THEN
196        WRITE(YRECFM,FMT='(A8,I1)') 'ROOTFRAC',JL
197      ELSE
198        WRITE(YRECFM,FMT='(A8,I2)') 'ROOTFRAC',JL          
199      ENDIF  
200      YCOMMENT='root fraction by layer'
201      YCOMMENTUNIT='-'
202      ZWORK(:,:)=XUNDEF
203      DO JJ=1,SIZE(XDG,1)
204         WHERE(JL<=NWG_LAYER(JJ,:).AND.NWG_LAYER(JJ,:)/=NUNDEF)
205               ZWORK(JJ,:)=XROOTFRAC(JJ,JL,:)
206         ENDWHERE
207      ENDDO
208      CALL WRITE_SURF_FIELD2D(HPROGRAM,ZWORK(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
209   END DO
210 !
211 !* SOC fraction for each layer
212 !
213   IF(CSOC=='SGH')THEN
214     DO JL=1,SIZE(XDG,2)
215      IF (JL<10) THEN
216        WRITE(YRECFM,FMT='(A7,I1)') 'FRACSOC',JL
217      ELSE
218        WRITE(YRECFM,FMT='(A7,I2)') 'FRACSOC',JL          
219      ENDIF  
220      YCOMMENT='SOC fraction by layer (-)'
221      CALL WRITE_SURF(HPROGRAM,YRECFM,XFRACSOC(:,JL),IRESP,HCOMMENT=YCOMMENT)
222     END DO
223   ENDIF
224 !
225 ENDIF        
226 !
227 !-------------------------------------------------------------------------------
228 !
229 DO JL=1,SIZE(XDG,2)
230    IF (JL<10) THEN
231      WRITE(YRECFM,FMT='(A4,I1)') 'WSAT',JL
232    ELSE
233      WRITE(YRECFM,FMT='(A4,I2)') 'WSAT',JL          
234    ENDIF  
235   YCOMMENT='soil porosity by layer (m3/m3)'
236   CALL WRITE_SURF(HPROGRAM,YRECFM,XWSAT(:,JL),IRESP,HCOMMENT=YCOMMENT)
237 ENDDO
238 !
239 DO JL=1,SIZE(XDG,2)
240    IF (JL<10) THEN
241      WRITE(YRECFM,FMT='(A3,I1)') 'WFC',JL
242    ELSE
243      WRITE(YRECFM,FMT='(A3,I2)') 'WFC',JL          
244    ENDIF  
245   YCOMMENT='field capacity by layer (m3/m3)'
246   CALL WRITE_SURF(HPROGRAM,YRECFM,XWFC(:,JL),IRESP,HCOMMENT=YCOMMENT)
247 ENDDO
248 !
249 DO JL=1,SIZE(XDG,2)
250    IF (JL<10) THEN
251      WRITE(YRECFM,FMT='(A5,I1)') 'WWILT',JL
252    ELSE
253      WRITE(YRECFM,FMT='(A5,I2)') 'WWILT',JL          
254    ENDIF  
255   YCOMMENT='wilting point by layer (m3/m3)'
256   CALL WRITE_SURF(HPROGRAM,YRECFM,XWWILT(:,JL),IRESP,HCOMMENT=YCOMMENT)
257 ENDDO     
258 !
259 !-------------------------------------------------------------------------------
260 ! For Earth System Model
261 IF(LFANOCOMPACT.AND..NOT.LPREP)THEN
262   CALL END_IO_SURF_n(HPROGRAM)
263   IF (LHOOK) CALL DR_HOOK('WRITE_DIAG_PGD_ISBA_N',1,ZHOOK_HANDLE)
264   RETURN
265 ENDIF
266 !
267 !-------------------------------------------------------------------------------
268 !
269 YRECFM='Z0REL'
270 YCOMMENT='orography roughness length (M)'
271 !
272  CALL WRITE_SURF(HPROGRAM,YRECFM,XZ0REL(:),IRESP,HCOMMENT=YCOMMENT)
273 !
274 !-------------------------------------------------------------------------------
275 !
276 !* Runoff soil ice depth for each patch
277 !
278 IF(CHORT=='SGH'.AND.CISBA/='DIF')THEN
279   YRECFM='DICE'
280   YCOMMENT='soil ice depth for runoff'
281   YCOMMENTUNIT='m'
282   CALL WRITE_SURF_FIELD2D(HPROGRAM,XD_ICE(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
283 ENDIF
284 !
285 !-------------------------------------------------------------------------------
286 !
287 !* Fraction of each vegetation type for each patch
288 !
289 DO JL=1,SIZE(XVEGTYPE_PATCH,2)
290   WRITE(YPAS,'(I2)') JL 
291   YLVLV=ADJUSTL(YPAS(:LEN_TRIM(YPAS)))
292   WRITE(YRECFM,FMT='(A9)') 'VEGTY_P'//YLVLV
293   YCOMMENT='fraction of each vegetation type for each patch'
294   YCOMMENTUNIT='-'
295   CALL WRITE_SURF_FIELD2D(HPROGRAM,XVEGTYPE_PATCH(:,JL,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
296 END DO
297 !-------------------------------------------------------------------------------
298 !
299 !* other surface parameters
300 !
301 YRECFM='RSMIN'
302 YCOMMENT='minimum stomatal resistance'
303 YCOMMENTUNIT='SM-1'
304 CALL WRITE_SURF_FIELD2D(HPROGRAM,XRSMIN(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
305 !
306 YRECFM='GAMMA'
307 YCOMMENT='coefficient for RSMIN calculation'
308 YCOMMENTUNIT='-'
309 CALL WRITE_SURF_FIELD2D(HPROGRAM,XGAMMA(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
310 !
311 YRECFM='CV'
312 YCOMMENT='vegetation thermal inertia coefficient'
313 YCOMMENTUNIT='-'
314 CALL WRITE_SURF_FIELD2D(HPROGRAM,XCV(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
315 !
316 YRECFM='RGL'
317 YCOMMENT='maximum solar radiation usable in photosynthesis'
318 YCOMMENTUNIT='-'
319 CALL WRITE_SURF_FIELD2D(HPROGRAM,XRGL(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
320 !
321 YRECFM='EMIS_ISBA'
322 YCOMMENT='surface emissivity'
323 YCOMMENTUNIT='-'
324 CALL WRITE_SURF_FIELD2D(HPROGRAM,XEMIS(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
325 !
326 YRECFM='WRMAX_CF'
327 YCOMMENT='coefficient for maximum water interception'
328 YCOMMENTUNIT='-'
329 CALL WRITE_SURF_FIELD2D(HPROGRAM,XWRMAX_CF(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
330 !
331 !-------------------------------------------------------------------------------
332 !
333 IF (LSURF_DIAG_ALBEDO) THEN
334 !
335 !* Soil albedos
336 !
337 !
338    YRECFM='ALBNIR_S'
339    YCOMMENT='soil near-infra-red albedo'
340    YCOMMENTUNIT='-'
341    CALL WRITE_SURF_FIELD2D(HPROGRAM,XALBNIR_SOIL(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
342 !
343 !-------------------------------------------------------------------------------
344 !
345    YRECFM='ALBVIS_S'
346    YCOMMENT='soil visible albedo'
347    YCOMMENTUNIT='-'
348    CALL WRITE_SURF_FIELD2D(HPROGRAM,XALBVIS_SOIL(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
349 !
350 !-------------------------------------------------------------------------------
351 !
352    YRECFM='ALBUV_S'
353    YCOMMENT='soil UV albedo'
354    YCOMMENTUNIT='-'
355    CALL WRITE_SURF_FIELD2D(HPROGRAM,XALBUV_SOIL(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
356 !
357 !-------------------------------------------------------------------------------
358 !
359 !* albedos
360 !
361    YRECFM='ALBNIR_ISBA'
362    YCOMMENT='total near-infra-red albedo'
363    YCOMMENTUNIT='-'
364    CALL WRITE_SURF_FIELD2D(HPROGRAM,XALBNIR(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
365 !
366 !-------------------------------------------------------------------------------
367 !
368    YRECFM='ALBVIS_ISBA'
369    YCOMMENT='total visible albedo'
370    YCOMMENTUNIT='-'
371    CALL WRITE_SURF_FIELD2D(HPROGRAM,XALBVIS(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
372 !
373 !-------------------------------------------------------------------------------
374 !
375    YRECFM='ALBUV_ISBA'
376    YCOMMENT='total UV albedo'
377    YCOMMENTUNIT='-'
378    CALL WRITE_SURF_FIELD2D(HPROGRAM,XALBUV(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
379 !
380 END IF
381 !
382 !-------------------------------------------------------------------------------
383 !
384 !* chemical soil resistances
385 !
386 IF (CCH_DRY_DEP=='WES89' .AND. NBEQ>0) THEN
387   YRECFM='SOILRC_SO2'
388   YCOMMENT='bare soil resistance for SO2'
389   YCOMMENTUNIT='?'
390   CALL WRITE_SURF_FIELD2D(HPROGRAM,XSOILRC_SO2(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
391   !
392   YRECFM='SOILRC_O3'
393   YCOMMENT='bare soil resistance for O3'
394   YCOMMENTUNIT='?'
395   CALL WRITE_SURF_FIELD2D(HPROGRAM,XSOILRC_O3(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
396 END IF
397 !
398 !-------------------------------------------------------------------------------
399 !
400 IF (LAGRIP .AND. (CPHOTO=='LAI' .OR. CPHOTO=='LST' .OR. CPHOTO=='NIT' .OR. CPHOTO=='NCB') ) THEN
401 !
402 !* seeding and reaping
403 !
404 !
405   YRECFM='TSEED'
406   YCOMMENT='date of seeding (-)'
407 !
408   CALL WRITE_SURF(HPROGRAM,YRECFM,TSEED(:,:),IRESP,HCOMMENT=YCOMMENT)
409 !
410   YRECFM='TREAP'
411   YCOMMENT='date of reaping (-)'
412 !
413   CALL WRITE_SURF(HPROGRAM,YRECFM,TREAP(:,:),IRESP,HCOMMENT=YCOMMENT)
414 !
415 !-------------------------------------------------------------------------------
416 !
417 !* irrigated fraction
418 !
419   YRECFM='IRRIG'
420   YCOMMENT='flag for irrigation (irrigation if >0.)'
421   YCOMMENTUNIT='-'
422 !
423   CALL WRITE_SURF_FIELD2D(HPROGRAM,XIRRIG(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
424 !
425 !-------------------------------------------------------------------------------
426 !
427 !* water supply for irrigation
428 !
429   YRECFM='WATSUP'
430   YCOMMENT='water supply during irrigation process'
431   YCOMMENTUNIT='mm'
432 !
433   CALL WRITE_SURF_FIELD2D(HPROGRAM,XWATSUP(:,:),YRECFM,YCOMMENT,YCOMMENTUNIT)
434 !
435 ENDIF
436 !-------------------------------------------------------------------------------
437 !         End of IO
438 !
439  CALL END_IO_SURF_n(HPROGRAM)
440 IF (LHOOK) CALL DR_HOOK('WRITE_DIAG_PGD_ISBA_N',1,ZHOOK_HANDLE)
441 !
442 !
443 END SUBROUTINE WRITE_DIAG_PGD_ISBA_n