Merge branch 'LIBTOOLS-master' into MNH-52X
[MNH-git_open_source-lfs.git] / src / SURFEX / read_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 READ_OCEAN_n(HPROGRAM)
7 !     #########################################
8 !
9 !!****  *READ_OCEAN_n* - read oceanic variables
10 !!
11 !!
12 !!    PURPOSE
13 !!    -------
14 !!
15 !!**  METHOD
16 !!    ------
17 !!
18 !!    EXTERNAL
19 !!    --------
20 !!
21 !!
22 !!    IMPLICIT ARGUMENTS
23 !!    ------------------
24 !!
25 !!    REFERENCE
26 !!    ---------
27 !!
28 !!
29 !!    AUTHOR
30 !!    ------
31 !!      C. Lebeaupin Brossier   *Meteo France*  
32 !!
33 !!    MODIFICATIONS
34 !!    -------------
35 !!      Original    04/2007 
36 !!      Modofied    07/2012, P. Le Moigne : CMO1D phasing
37 !-------------------------------------------------------------------------------
38 !
39 !*       0.    DECLARATIONS
40 !              ------------
41 !
42 USE MODD_SURF_PAR, ONLY : XUNDEF
43 USE MODD_OCEAN_n
44 USE MODD_DIAG_OCEAN_n
45 USE MODD_OCEAN_CSTS, ONLY : NOCKMIN,NOCKMAX
46 !
47 !
48 USE MODI_READ_SURF
49 USE MODI_OCEAN_MERCATORVERGRID
50 !
51 USE MODD_OCEAN_REL_n
52 !
53 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
54 USE PARKIND1  ,ONLY : JPRB
55 !
56 USE MODI_GET_TYPE_DIM_n
57 !
58 IMPLICIT NONE
59 !
60 !*       0.1   Declarations of arguments
61 !              -------------------------
62 !
63  CHARACTER(LEN=6),  INTENT(IN)  :: HPROGRAM ! calling program
64 !
65 !*       0.2   Declarations of local variables
66 !              -------------------------------
67 !
68 INTEGER           :: ILU          ! 1D physical dimension
69 !
70 INTEGER           :: IRESP          ! Error code after redding
71 !
72  CHARACTER(LEN=4)  :: YLVL
73 !
74  CHARACTER(LEN=12) :: YRECFM         ! Name of the article to be read
75  CHARACTER(LEN=14) :: YFORM          ! Writing format
76 REAL, DIMENSION(:),ALLOCATABLE  :: ZWORK      ! 1D array to write data in file
77 !
78 INTEGER :: JLEVEL ! loop counter on oceanic levels
79 INTEGER :: J      ! loop counter on sea grid points
80 !
81 INTEGER           :: IVERSION       ! surface version
82 REAL(KIND=JPRB) :: ZHOOK_HANDLE
83 !
84 !-------------------------------------------------------------------------------
85 IF (LHOOK) CALL DR_HOOK('READ_OCEAN_N',0,ZHOOK_HANDLE)
86 NOCTCOUNT=0
87 !
88 YRECFM='VERSION'
89  CALL READ_SURF(HPROGRAM,YRECFM,IVERSION,IRESP)
90 !
91 !* flag to use or not Ocean model
92 !
93 IF (IVERSION<=3) THEN
94    LMERCATOR=.FALSE.
95 ELSE
96    YRECFM='SEA_OCEAN'
97    CALL READ_SURF(HPROGRAM,YRECFM,LMERCATOR,IRESP)
98 ENDIF
99 !
100 IF (.NOT. LMERCATOR) THEN
101   ALLOCATE(XSEAT(0,0))
102   ALLOCATE(XSEAT_REL(0,0))
103   ALLOCATE(XSEAS(0,0))
104   ALLOCATE(XSEAS_REL(0,0))
105   ALLOCATE(XSEAU_REL(0,0))
106   ALLOCATE(XSEAV_REL(0,0))
107   ALLOCATE(XSEAU(0,0))
108   ALLOCATE(XSEAV(0,0))
109   ALLOCATE(XSEAE(0,0))
110   ALLOCATE(XSEABATH(0,0))
111   ALLOCATE(XSEAHMO(0))
112   ALLOCATE(XLE        (0,0))
113   ALLOCATE(XLK        (0,0))
114   ALLOCATE(XKMEL      (0,0))
115   ALLOCATE(XKMELM     (0,0))
116   ALLOCATE(XSEATEND   (0))
117   ALLOCATE(XDTFSOL(0,0))
118   ALLOCATE(XDTFNSOL(0))
119   IF (LHOOK) CALL DR_HOOK('READ_OCEAN_N',1,ZHOOK_HANDLE)
120   RETURN
121 ENDIF
122 !
123 !-------------------------------------------------------------------------------
124 !$OMP SINGLE
125  CALL OCEAN_MERCATORVERGRID
126 !$OMP END SINGLE
127 !
128 ! Relaxation time and logical
129 YRECFM='TAU_REL_OC'
130  CALL READ_SURF(HPROGRAM,YRECFM,XTAU_REL,IRESP)
131 !
132 YRECFM='LREL_CUR_OC'
133  CALL READ_SURF(HPROGRAM,YRECFM,LREL_CUR,IRESP)
134
135 YRECFM='LREL_TS_OC'
136  CALL READ_SURF(HPROGRAM,YRECFM,LREL_TS,IRESP)
137 YRECFM='LFLX_NULL_OC'
138  CALL READ_SURF(HPROGRAM,YRECFM,LFLUX_NULL,IRESP)
139 YRECFM='LFLX_CORR_OC'
140  CALL READ_SURF(HPROGRAM,YRECFM,LFLX_CORR,IRESP)
141 YRECFM='CORR_FLX_OC'
142  CALL READ_SURF(HPROGRAM,YRECFM,XQCORR,IRESP)
143 YRECFM='LDIAPYC_OC'
144  CALL READ_SURF(HPROGRAM,YRECFM,LDIAPYCNAL,IRESP)
145 !
146 !* 1D physical dimension
147 !
148 YRECFM='SIZE_SEA'
149  CALL GET_TYPE_DIM_n('SEA   ',ILU)
150 !
151 !*       2.     Prognostic fields:
152 !               -----------------
153 !
154 ALLOCATE(ZWORK(ILU))
155 !* oceanic temperature
156 !
157 ALLOCATE(XSEAT(ILU,NOCKMIN:NOCKMAX))
158 !
159 DO JLEVEL=NOCKMIN+1,NOCKMAX
160   WRITE(YLVL,'(I4)') JLEVEL
161   YRECFM='TEMP_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
162   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
163   XSEAT(:,JLEVEL)=ZWORK(:)
164 END DO
165 XSEAT(:,NOCKMIN)=XSEAT(:,NOCKMIN+1)
166 !
167 !* relaxation profile for oceanic temperature
168 !
169 ALLOCATE(XSEAT_REL(ILU,NOCKMIN:NOCKMAX))
170 !
171 DO JLEVEL=NOCKMIN+1,NOCKMAX
172   WRITE(YLVL,'(I4)') JLEVEL
173   YRECFM='T_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
174   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
175   XSEAT_REL(:,JLEVEL)=ZWORK(:)
176 END DO
177 XSEAT_REL(:,NOCKMIN)=XSEAT_REL(:,NOCKMIN+1)
178 !
179 !* oceanic salinity
180 !
181 ALLOCATE(XSEAS(ILU,NOCKMIN:NOCKMAX))
182 !
183 DO JLEVEL=NOCKMIN+1,NOCKMAX
184   WRITE(YLVL,'(I4)') JLEVEL
185   YRECFM='SALT_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
186   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
187   XSEAS(:,JLEVEL)=ZWORK(:)
188 END DO
189 XSEAS(:,NOCKMIN)=XSEAS(:,NOCKMIN+1)
190 !
191 !* oceanic salinity profile of relaxation
192 !
193 ALLOCATE(XSEAS_REL(ILU,NOCKMIN:NOCKMAX))
194 !
195 DO JLEVEL=NOCKMIN+1,NOCKMAX
196   WRITE(YLVL,'(I4)') JLEVEL
197   YRECFM='S_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
198   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
199   XSEAS_REL(:,JLEVEL)=ZWORK(:)
200 END DO
201 XSEAS_REL(:,NOCKMIN)=XSEAS_REL(:,NOCKMIN+1)
202 !
203 !* oceanic current
204 !
205 ALLOCATE(XSEAU_REL(ILU,NOCKMIN:NOCKMAX))
206 ALLOCATE(XSEAV_REL(ILU,NOCKMIN:NOCKMAX))
207 !
208 DO JLEVEL=NOCKMIN+1,NOCKMAX
209   WRITE(YLVL,'(I4)') JLEVEL
210   YRECFM='U_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
211   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
212   XSEAU_REL(:,JLEVEL)=ZWORK(:)
213 END DO
214 XSEAU_REL(:,NOCKMIN)=XSEAU_REL(:,NOCKMIN+1)
215 !
216 DO JLEVEL=NOCKMIN+1,NOCKMAX
217   WRITE(YLVL,'(I4)') JLEVEL
218   YRECFM='V_OC_REL'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
219   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
220   XSEAV_REL(:,JLEVEL)=ZWORK(:)
221 END DO
222 XSEAV_REL(:,NOCKMIN)=XSEAV_REL(:,NOCKMIN+1)
223 !
224 ALLOCATE(XSEAU(ILU,NOCKMIN:NOCKMAX))
225 ALLOCATE(XSEAV(ILU,NOCKMIN:NOCKMAX))
226 !
227 DO JLEVEL=NOCKMIN+1,NOCKMAX
228   WRITE(YLVL,'(I4)') JLEVEL
229   YRECFM='UCUR_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
230   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
231   XSEAU(:,JLEVEL)=ZWORK(:)
232 END DO
233 DO JLEVEL=NOCKMIN+1,NOCKMAX
234   WRITE(YLVL,'(I4)') JLEVEL
235   YRECFM='VCUR_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
236   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
237   XSEAV(:,JLEVEL)=ZWORK(:)
238 END DO
239 XSEAU(:,NOCKMIN)=XSEAU(:,NOCKMIN+1)
240 XSEAV(:,NOCKMIN)=XSEAV(:,NOCKMIN+1)
241 !
242 !* oceanic turbulent kinetic energy
243 !
244 ALLOCATE(XSEAE(ILU,NOCKMIN:NOCKMAX))
245 !
246 DO JLEVEL=NOCKMIN+1,NOCKMAX
247   WRITE(YLVL,'(I4)') JLEVEL
248   YRECFM='TKE_OC'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
249   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
250   XSEAE(:,JLEVEL)=ZWORK(:)
251 END DO
252 XSEAE(:,NOCKMIN)=XSEAE(:,NOCKMIN+1)
253 !
254 !
255 !-------------------------------------------------------------------------------
256 !
257 !*       4.     Semi-prognostic fields:
258 !               ----------------------
259 !
260 !* bathymetry indice
261 !
262 ALLOCATE(XSEABATH(ILU,NOCKMIN:NOCKMAX))
263 !
264 DO JLEVEL=NOCKMIN+1,NOCKMAX
265   WRITE(YLVL,'(I4)') JLEVEL
266   YRECFM='SEAINDBATH'//ADJUSTL(YLVL(:LEN_TRIM(YLVL)))
267   CALL READ_SURF(HPROGRAM,YRECFM,ZWORK(:),IRESP)
268   XSEABATH(:,JLEVEL)=ZWORK(:)
269 END DO
270 XSEABATH(:,NOCKMIN)=1.
271 !
272 !-------------------------------------------------------------------------------
273 !Complete undefined values for the oceanic 1D model convergence
274 DO J=1,ILU
275   DO JLEVEL=NOCKMIN+2,NOCKMAX
276     IF (XSEABATH(J,JLEVEL)==0.) THEN
277       XSEAT(J,JLEVEL)=XSEAT(J,JLEVEL-1)
278       XSEAS(J,JLEVEL)=XSEAS(J,JLEVEL-1)
279       XSEAU(J,JLEVEL)=XSEAU(J,JLEVEL-1)
280       XSEAV(J,JLEVEL)=XSEAV(J,JLEVEL-1)
281       XSEAE(J,JLEVEL)=XSEAE(J,JLEVEL-1)
282       !
283       XSEAT_REL(J,JLEVEL)=XSEAT_REL(J,JLEVEL-1)
284       XSEAS_REL(J,JLEVEL)=XSEAS_REL(J,JLEVEL-1)
285       XSEAU_REL(J,JLEVEL)=XSEAU_REL(J,JLEVEL-1)
286       XSEAV_REL(J,JLEVEL)=XSEAV_REL(J,JLEVEL-1)      
287     ENDIF
288   ENDDO
289 ENDDO
290 !
291 DEALLOCATE(ZWORK)
292 !-------------------------------------------------------------------------------
293 ALLOCATE(XSEAHMO(ILU))
294 YRECFM='SEA_HMO'
295  CALL READ_SURF(HPROGRAM,YRECFM,XSEAHMO(:),IRESP)
296 !
297 !-------------------------------------------------------------------------------
298 ALLOCATE(XLE        (SIZE(XSEAT,1),NOCKMIN:NOCKMAX))
299 ALLOCATE(XLK        (SIZE(XSEAT,1),NOCKMIN:NOCKMAX))
300 ALLOCATE(XKMEL      (SIZE(XSEAT,1),NOCKMIN:NOCKMAX))
301 ALLOCATE(XKMELM     (SIZE(XSEAT,1),NOCKMIN:NOCKMAX))
302 XLE(:,:)    =XUNDEF
303 XLK(:,:)    =XUNDEF
304 XKMEL(:,:)  =XUNDEF
305 XKMELM(:,:) =XUNDEF
306 !
307 ALLOCATE(XSEATEND   (SIZE(XSEAT,1)))
308 XSEATEND(:) =XUNDEF
309 !
310 ALLOCATE(XDTFSOL(ILU,NOCKMIN:NOCKMAX))
311 ALLOCATE(XDTFNSOL(ILU))
312 !
313 XDTFSOL(:,:) = XUNDEF 
314 XDTFNSOL(:) = XUNDEF 
315 !
316 IF (LHOOK) CALL DR_HOOK('READ_OCEAN_N',1,ZHOOK_HANDLE)
317 !
318 !------------------------------------------------------------------------------
319 END SUBROUTINE READ_OCEAN_n