Juan 8/12/2016: add management of LEN_HREC in MNH & SURFEX
[MNH-git_open_source-lfs.git] / src / SURFEX / writesurf_oceann.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 WRITESURF_OCEAN_n(HPROGRAM)
7 !     ########################################
8 !
9 !!****  *WRITE_OCEAN_n* - writes OCEAN 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 !!      C. Lebeaupin Brossier   *Meteo France*  
31 !!
32 !!    MODIFICATIONS
33 !!    -------------
34 !!      Original    03/2007
35 !!      Modified    07/2012, P. Le Moigne : CMO1D phasing
36 !-------------------------------------------------------------------------------
37 !
38 !*       0.    DECLARATIONS
39 !              ------------
40 !
41 USE MODD_OCEAN_n, ONLY : XSEAT, XSEAS, XSEAU, XSEAV, XSEAE, XSEABATH,&
42                          XSEAHMO, LMERCATOR, XDTFSOL, XDTFNSOL, XKMEL
43 USE MODD_OCEAN_GRID_n
44 !
45 USE MODD_OCEAN_CSTS
46 !
47 USE MODI_WRITE_SURF
48 !
49 USE MODD_OCEAN_REL_n, ONLY : XSEAT_REL, XSEAS_REL,XTAU_REL,     & 
50                            LREL_CUR,LREL_TS,LFLUX_NULL,XQCORR,LFLX_CORR,&
51                            XSEAU_REL,XSEAV_REL,LDIAPYCNAL
52 !
53 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
54 USE PARKIND1  ,ONLY : JPRB
55 !
56 IMPLICIT NONE
57 !
58 !*       0.1   Declarations of arguments
59 !              -------------------------
60 !
61  CHARACTER(LEN=6),  INTENT(IN)  :: HPROGRAM ! program calling
62
63 !
64 !*       0.2   Declarations of local variables
65 !              -------------------------------
66 !
67 INTEGER           :: IRESP          ! IRESP  : return-code if a problem appears
68  CHARACTER(LEN=LEN_HREC) :: YRECFM         ! Name of the article to be read
69  CHARACTER(LEN=4 ) :: YLVL
70  CHARACTER(LEN=100):: YCOMMENT       ! Comment string
71  CHARACTER(LEN=14) :: YFORM          ! Writing format
72 !
73 INTEGER :: JLEVEL ! loop counter on oceanic levels
74 !
75 REAL, DIMENSION(SIZE(XSEAT,1)) :: ZWORK ! 1D array to write data in file
76 REAL(KIND=JPRB) :: ZHOOK_HANDLE
77 !
78 !-------------------------------------------------------------------------------
79 IF (LHOOK) CALL DR_HOOK('WRITESURF_OCEAN_N',0,ZHOOK_HANDLE)
80 !
81 !* flag to define if OCEAN is used
82 !
83 YRECFM='SEA_OCEAN'
84 YCOMMENT='flag to use OCEAN model'
85  CALL WRITE_SURF(HPROGRAM,YRECFM,LMERCATOR,IRESP,HCOMMENT=YCOMMENT)
86 !
87 IF (.NOT. LMERCATOR .AND. LHOOK) CALL DR_HOOK('WRITESURF_OCEAN_N',1,ZHOOK_HANDLE)
88 IF (.NOT. LMERCATOR) RETURN
89 !
90 ! Relaxation time  
91 YRECFM='TAU_REL_OC'
92 YCOMMENT='Relaxation time of ocean model (s)'
93  CALL WRITE_SURF(HPROGRAM,YRECFM,XTAU_REL,IRESP,HCOMMENT=YCOMMENT)
94 !
95 YRECFM='LREL_CUR_OC'
96 YCOMMENT='FLAG FOR RELAXATION ON CURRENT'
97  CALL WRITE_SURF(HPROGRAM,YRECFM,LREL_CUR,IRESP,HCOMMENT=YCOMMENT)
98 !
99 YRECFM='LREL_TS_OC'
100 YCOMMENT='FLAG FOR RELAXATION ON T,S'
101  CALL WRITE_SURF(HPROGRAM,YRECFM,LREL_TS,IRESP,HCOMMENT=YCOMMENT)
102 !
103 YRECFM='LFLX_NULL_OC'
104 YCOMMENT='FLAG FOR ZERO FLUX '
105  CALL WRITE_SURF(HPROGRAM,YRECFM,LFLUX_NULL,IRESP,HCOMMENT=YCOMMENT)
106
107 YRECFM='LFLX_CORR_OC'
108 YCOMMENT='FLAG FOR FLUX CORRECTION '
109  CALL WRITE_SURF(HPROGRAM,YRECFM,LFLX_CORR,IRESP,HCOMMENT=YCOMMENT)
110
111 YRECFM='CORR_FLX_OC'
112 YCOMMENT='FLUX CORRECTION COEFF (W/M2/K)'
113  CALL WRITE_SURF(HPROGRAM,YRECFM,XQCORR,IRESP,HCOMMENT=YCOMMENT)
114 !
115 YRECFM='LDIAPYC_OC'
116 YCOMMENT='FLAG FOR DIAPYCNAL MIXING '
117  CALL WRITE_SURF(HPROGRAM,YRECFM,LDIAPYCNAL,IRESP,HCOMMENT=YCOMMENT)
118 !
119 !-------------------------------------------------------------------------------
120 !
121 !*       3.     Prognostic fields:
122 !               -----------------
123 !
124 ! tendance surface liee au flux non solaire
125 !
126   JLEVEL=1
127   WRITE(YLVL,'(I4)') JLEVEL
128   YRECFM='DTFNSOL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
129   YFORM='(A11,I1.1,A5)'
130   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
131   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_DTFNSOL',JLEVEL,' (°C)'
132   ZWORK=XDTFNSOL(:)
133   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)
134 !
135 !* oceanic temperature
136 !
137 DO JLEVEL=NOCKMIN+1,NOCKMAX
138   WRITE(YLVL,'(I4)') JLEVEL
139   YRECFM='TEMP_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
140   YFORM='(A11,I1.1,A5)'
141   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
142   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_TEMP_OC',JLEVEL,' (°C)'
143   ZWORK=XSEAT(:,JLEVEL)
144   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)
145 END DO
146 !
147 ! Budget terms lie au flux solaire
148 DO JLEVEL=NOCKMIN+1,NOCKMAX
149   WRITE(YLVL,'(I4)') JLEVEL
150   YRECFM='DTFSOL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
151   YFORM='(A11,I1.1,A5)'
152   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
153   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_DTFSOL',JLEVEL,' (°C/s)'
154   ZWORK=XDTFSOL(:,JLEVEL)
155   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)
156 END DO
157
158 !
159 !* relaxation profile for oceanic temperature
160 DO JLEVEL=NOCKMIN+1,NOCKMAX
161   WRITE(YLVL,'(I4)') JLEVEL
162   YRECFM='T_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
163   YFORM='(A11,I1.1,A5)'
164   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
165   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_T_OC_REL',JLEVEL,' (°C)'
166   ZWORK=XSEAT_REL(:,JLEVEL)
167   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)  
168 END DO
169 !
170 !* oceanic salinity
171 !
172 DO JLEVEL=NOCKMIN+1,NOCKMAX
173   WRITE(YLVL,'(I4)') JLEVEL
174   YRECFM='SALT_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
175   YFORM='(A11,I1.1,A5)'
176   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
177   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_SALT_OC',JLEVEL,'(psu)'
178   ZWORK=XSEAS(:,JLEVEL)
179   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)
180 END DO
181 !
182 !* oceanic salinity profile of relxation
183 !
184 DO JLEVEL=NOCKMIN+1,NOCKMAX
185   WRITE(YLVL,'(I4)') JLEVEL
186   YRECFM='S_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
187   YFORM='(A11,I1.1,A5)'
188   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
189   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_S_OC_REL',JLEVEL,'(psu)'
190   ZWORK=XSEAS_REL(:,JLEVEL)
191   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)  
192 END DO
193 !
194 !* oceanic zonal current
195 !
196 DO JLEVEL=NOCKMIN+1,NOCKMAX
197   WRITE(YLVL,'(I4)') JLEVEL
198   YRECFM='U_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
199   YFORM='(A11,I1.1,A5)'
200   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
201   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_U_OC_REL',JLEVEL,' M/S'
202   ZWORK=XSEAU_REL(:,JLEVEL)
203   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)
204 END DO
205 !
206 DO JLEVEL=NOCKMIN+1,NOCKMAX
207   WRITE(YLVL,'(I4)') JLEVEL
208   YRECFM='UCUR_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
209   YFORM='(A11,I1.1,A5)'
210   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
211   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_UCUR_OC',JLEVEL,' (m/s)'
212   ZWORK=XSEAU(:,JLEVEL)
213   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)  
214 END DO
215 !
216 !* oceanic meridian current
217 !
218 DO JLEVEL=NOCKMIN+1,NOCKMAX
219   WRITE(YLVL,'(I4)') JLEVEL
220   YRECFM='V_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
221   YFORM='(A11,I1.1,A5)'
222   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
223   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_V_OC_REL',JLEVEL,' M/S'
224   ZWORK=XSEAV_REL(:,JLEVEL)
225   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)
226 END DO
227 !
228 DO JLEVEL=NOCKMIN+1,NOCKMAX
229   WRITE(YLVL,'(I4)') JLEVEL  
230   YRECFM='VCUR_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
231   YFORM='(A11,I1.1,A5)'
232   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
233   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_VCUR_OC',JLEVEL,'(m/s)'
234   ZWORK=XSEAV(:,JLEVEL)
235   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)
236 END DO
237 !
238 !* oceanic turbulent kinetic energy
239 !
240 DO JLEVEL=NOCKMIN+1,NOCKMAX
241   WRITE(YLVL,'(I4)') JLEVEL
242   YRECFM='TKE_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
243   YFORM='(A11,I1.1,A5)'
244   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
245   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_TKE_OC',JLEVEL,' (J)'
246   ZWORK=XSEAE(:,JLEVEL)
247   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)  
248 END DO
249 !
250 !* oceanic mixing coefficient
251 !
252 DO JLEVEL=NOCKMIN+1,NOCKMAX
253   WRITE(YLVL,'(I4)') JLEVEL
254   YRECFM='KMEL_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
255   YFORM='(A11,I1.1,A5)'
256   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
257   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_KMEL_OC',JLEVEL,' (m2/s2)'
258   ZWORK=XKMEL(:,JLEVEL)
259   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)  
260 END DO
261
262 !
263 !-------------------------------------------------------------------------------
264 !
265 !*       4.     Semi-prognostic fields:
266 !               ----------------------
267 !* bathymetry indice
268 !
269 DO JLEVEL=NOCKMIN+1,NOCKMAX
270   WRITE(YLVL,'(I4)') JLEVEL
271   YRECFM='SEAINDBATH'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
272   YFORM='(A11,I1.1,A5)'
273   IF (JLEVEL >= 10)  YFORM='(A11,I2.2,A5)'
274   WRITE(YCOMMENT,FMT=YFORM) 'X_Y_SEAINDBATH',JLEVEL,' (J)'
275   ZWORK=XSEABATH(:,JLEVEL)
276   CALL WRITE_SURF(HPROGRAM,YRECFM,ZWORK,IRESP,HCOMMENT=YCOMMENT)  
277 END DO
278 !
279 !-------------------------------------------------------------------------------
280 !Sea Surface Salinity
281 !
282 YRECFM='SSS'
283 YCOMMENT='SSS'
284  CALL WRITE_SURF(HPROGRAM,YRECFM,XSEAS(:,NOCKMIN),IRESP,HCOMMENT=YCOMMENT)
285 !-------------------------------------------------------------------------------
286 !Sea Surface Salinity
287 !
288 YRECFM='SEA_HMO'
289 YCOMMENT='X_Y_'
290  CALL WRITE_SURF(HPROGRAM,YRECFM,XSEAHMO(:),IRESP,HCOMMENT=YCOMMENT)
291 !-------------------------------------------------------------------------------
292 IF (LHOOK) CALL DR_HOOK('WRITESURF_OCEAN_N',1,ZHOOK_HANDLE)
293 !
294 END SUBROUTINE WRITESURF_OCEAN_n