source: branches/UKMO/r6232_tracer_advection/NEMOGCM/TOOLS/NESTING/src/agrif_readwrite.f90 @ 9295

Last change on this file since 9295 was 9295, checked in by jcastill, 3 years ago

Remove svn keywords

File size: 26.4 KB
Line 
1!************************************************************************
2! Fortran 95 OPA Nesting tools                  *
3!                          *
4!     Copyright (C) 2005 Florian Lemarié (Florian.Lemarie@imag.fr)   *
5!                          *
6!************************************************************************
7!
8MODULE agrif_readwrite
9  !
10  USE agrif_types
11  !
12  IMPLICIT NONE
13  !
14CONTAINS
15  !
16  !************************************************************************
17  !                           *
18  ! MODULE  AGRIF_READWRITE                  *
19  !                           *
20  ! module containing subroutine used for :           *
21  !   - Coordinates files reading/writing          *
22  !   - Bathymetry files reading/writing (meter and levels)    *
23  !   - Naming of child grid files              *
24  !                           *
25  !************************************************************************
26  !       
27  !*****************************************************
28  !   function Read_Coordinates(name,Grid)
29  !*****************************************************
30
31  INTEGER FUNCTION Read_Coordinates(name,Grid,Pacifique)
32    !
33    USE io_netcdf
34    !   
35    !  file name to open
36    !
37    CHARACTER(*) name
38    LOGICAL,OPTIONAL :: Pacifique
39    !
40    TYPE(Coordinates) :: Grid
41    !     
42    CALL Read_Ncdf_var('glamt',name,Grid%glamt)
43    CALL Read_Ncdf_var('glamu',name,Grid%glamu)
44    CALL Read_Ncdf_var('glamv',name,Grid%glamv)
45    CALL Read_Ncdf_var('glamf',name,Grid%glamf)
46    CALL Read_Ncdf_var('gphit',name,Grid%gphit)
47    CALL Read_Ncdf_var('gphiu',name,Grid%gphiu)
48    CALL Read_Ncdf_var('gphiv',name,Grid%gphiv)
49    CALL Read_Ncdf_var('gphif',name,Grid%gphif)
50    CALL Read_Ncdf_var('e1t',name,Grid%e1t)
51    CALL Read_Ncdf_var('e1u',name,Grid%e1u)
52    CALL Read_Ncdf_var('e1v',name,Grid%e1v)
53    CALL Read_Ncdf_var('e1f',name,Grid%e1f)
54    CALL Read_Ncdf_var('e2t',name,Grid%e2t)
55    CALL Read_Ncdf_var('e2u',name,Grid%e2u)
56    CALL Read_Ncdf_var('e2v',name,Grid%e2v)
57    CALL Read_Ncdf_var('e2f',name,Grid%e2f)
58    CALL Read_Ncdf_var('nav_lon',name,Grid%nav_lon)
59    CALL Read_Ncdf_var('nav_lat',name,Grid%nav_lat)       
60    !
61    IF( PRESENT(Pacifique) )THEN
62       IF ( Grid%glamt(1,1) > Grid%glamt(nxfin,nyfin) ) THEN           
63       Pacifique = .TRUE.
64       WHERE ( Grid%glamt < 0 )
65          Grid%glamt = Grid%glamt + 360.
66       END WHERE
67       WHERE ( Grid%glamf < 0 )
68          Grid%glamf = Grid%glamf + 360.
69       END WHERE
70       WHERE ( Grid%glamu < 0 )
71          Grid%glamu = Grid%glamu + 360.
72       END WHERE
73       WHERE ( Grid%glamv < 0 )
74          Grid%glamv = Grid%glamv + 360.
75       END WHERE
76       WHERE ( Grid%nav_lon < 0 )
77          Grid%nav_lon = Grid%nav_lon + 360.
78       END WHERE
79       ENDIF
80    ENDIF
81    !           
82    WRITE(*,*) ' '
83    WRITE(*,*) 'Reading coordinates file: ',name
84    WRITE(*,*) ' '
85    !     
86    Read_Coordinates = 1
87    !     
88  END FUNCTION Read_Coordinates
89
90  !*****************************************************
91  !   function Read_Coordinates(name,Grid)
92  !*****************************************************
93
94  INTEGER FUNCTION Read_Local_Coordinates(name,Grid,strt,cnt)
95    !
96    USE io_netcdf
97    !   
98    !  file name to open
99    !
100    CHARACTER(*) name
101    INTEGER, DIMENSION(2) :: strt,cnt
102    !
103    TYPE(Coordinates) :: Grid
104    !     
105    CALL Read_Ncdf_var('glamt',name,Grid%glamt,strt,cnt)
106    CALL Read_Ncdf_var('glamu',name,Grid%glamu,strt,cnt)
107    CALL Read_Ncdf_var('glamv',name,Grid%glamv,strt,cnt)
108    CALL Read_Ncdf_var('glamf',name,Grid%glamf,strt,cnt)
109    CALL Read_Ncdf_var('gphit',name,Grid%gphit,strt,cnt)
110    CALL Read_Ncdf_var('gphiu',name,Grid%gphiu,strt,cnt)
111    CALL Read_Ncdf_var('gphiv',name,Grid%gphiv,strt,cnt)
112    CALL Read_Ncdf_var('gphif',name,Grid%gphif,strt,cnt)
113    CALL Read_Ncdf_var('e1t',name,Grid%e1t,strt,cnt)
114    CALL Read_Ncdf_var('e1u',name,Grid%e1u,strt,cnt)
115    CALL Read_Ncdf_var('e1v',name,Grid%e1v,strt,cnt)
116    CALL Read_Ncdf_var('e1f',name,Grid%e1f,strt,cnt)
117    CALL Read_Ncdf_var('e2t',name,Grid%e2t,strt,cnt)
118    CALL Read_Ncdf_var('e2u',name,Grid%e2u,strt,cnt)
119    CALL Read_Ncdf_var('e2v',name,Grid%e2v,strt,cnt)
120    CALL Read_Ncdf_var('e2f',name,Grid%e2f,strt,cnt)
121    CALL Read_Ncdf_var('nav_lon',name,Grid%nav_lon,strt,cnt)
122    CALL Read_Ncdf_var('nav_lat',name,Grid%nav_lat,strt,cnt)       
123    !
124    WRITE(*,*) ' '
125    WRITE(*,*) 'Reading coordinates file: ',name
126    WRITE(*,*) ' '
127    !     
128    Read_Local_Coordinates = 1
129    !     
130  END FUNCTION Read_Local_Coordinates
131
132  !*****************************************************
133  !   function Write_Coordinates(name,Grid)
134  !*****************************************************
135
136  INTEGER FUNCTION Write_Coordinates(name,Grid)
137    !
138    USE io_netcdf
139    CHARACTER(*) name
140    TYPE(Coordinates) :: Grid
141    INTEGER :: status,ncid
142    CHARACTER(len=1),DIMENSION(2) :: dimnames
143    !
144    status = nf90_create(name,NF90_WRITE,ncid)
145    status = nf90_close(ncid)
146    !           
147    dimnames = (/ 'x','y' /)
148    CALL Write_Ncdf_dim(dimnames(1),name,nxfin)
149    CALL Write_Ncdf_dim(dimnames(2),name,nyfin)
150    !     
151    CALL Write_Ncdf_var('nav_lon',dimnames,name,Grid%nav_lon,'float')     
152    CALL Write_Ncdf_var('nav_lat',dimnames,name,Grid%nav_lat,'float')
153    !
154    CALL Write_Ncdf_var('glamt',dimnames,name,Grid%glamt,'double')
155    CALL Write_Ncdf_var('glamu',dimnames,name,Grid%glamu,'double')
156    CALL Write_Ncdf_var('glamv',dimnames,name,Grid%glamv,'double')
157    CALL Write_Ncdf_var('glamf',dimnames,name,Grid%glamf,'double')
158    CALL Write_Ncdf_var('gphit',dimnames,name,Grid%gphit,'double')
159    CALL Write_Ncdf_var('gphiu',dimnames,name,Grid%gphiu,'double')
160    CALL Write_Ncdf_var('gphiv',dimnames,name,Grid%gphiv,'double')
161    CALL Write_Ncdf_var('gphif',dimnames,name,Grid%gphif,'double')     
162    CALL Write_Ncdf_var('e1t',dimnames,name,Grid%e1t,'double')     
163    CALL Write_Ncdf_var('e1u',dimnames,name,Grid%e1u,'double')     
164    CALL Write_Ncdf_var('e1v',dimnames,name,Grid%e1v,'double')     
165    CALL Write_Ncdf_var('e1f',dimnames,name,Grid%e1f,'double')
166    CALL Write_Ncdf_var('e2t',dimnames,name,Grid%e2t,'double')
167    CALL Write_Ncdf_var('e2u',dimnames,name,Grid%e2u,'double')
168    CALL Write_Ncdf_var('e2v',dimnames,name,Grid%e2v,'double')
169    CALL Write_Ncdf_var('e2f',dimnames,name,Grid%e2f,'double')
170    !     
171    CALL Copy_Ncdf_att('nav_lon',TRIM(parent_coordinate_file),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon))
172    CALL Copy_Ncdf_att('nav_lat',TRIM(parent_coordinate_file),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat))
173    CALL Copy_Ncdf_att('glamt',TRIM(parent_coordinate_file),name)
174    CALL Copy_Ncdf_att('glamu',TRIM(parent_coordinate_file),name)
175    CALL Copy_Ncdf_att('glamv',TRIM(parent_coordinate_file),name)
176    CALL Copy_Ncdf_att('glamf',TRIM(parent_coordinate_file),name)
177    CALL Copy_Ncdf_att('gphit',TRIM(parent_coordinate_file),name)
178    CALL Copy_Ncdf_att('gphiu',TRIM(parent_coordinate_file),name)
179    CALL Copy_Ncdf_att('gphiv',TRIM(parent_coordinate_file),name)
180    CALL Copy_Ncdf_att('gphif',TRIM(parent_coordinate_file),name)
181    CALL Copy_Ncdf_att('e1t',TRIM(parent_coordinate_file),name)
182    CALL Copy_Ncdf_att('e1u',TRIM(parent_coordinate_file),name)
183    CALL Copy_Ncdf_att('e1v',TRIM(parent_coordinate_file),name)
184    CALL Copy_Ncdf_att('e1f',TRIM(parent_coordinate_file),name)
185    CALL Copy_Ncdf_att('e2t',TRIM(parent_coordinate_file),name)
186    CALL Copy_Ncdf_att('e2u',TRIM(parent_coordinate_file),name)
187    CALL Copy_Ncdf_att('e2v',TRIM(parent_coordinate_file),name)
188    CALL Copy_Ncdf_att('e2f',TRIM(parent_coordinate_file),name)           
189    !
190    WRITE(*,*) ' '
191    WRITE(*,*) 'Writing coordinates file: ',name
192    WRITE(*,*) ' '
193    !
194    Write_Coordinates = 1
195    !     
196  END FUNCTION Write_Coordinates
197  !
198  !
199  !
200  !*****************************************************
201  !   function Read_Bathy_level(name,Grid)
202  !*****************************************************
203  !
204  INTEGER FUNCTION Read_Bathy_level(name,Grid)
205    !
206    USE io_netcdf
207    !     
208    CHARACTER(*) name
209    TYPE(Coordinates) :: Grid
210    !
211    CALL Read_Ncdf_var('mbathy',name,Grid%Bathy_level)   
212    !
213    WRITE(*,*) ' '
214    WRITE(*,*) 'Reading bathymetry file: ',name
215    WRITE(*,*) ' '     
216    !
217    Read_Bathy_level = 1
218    !     
219  END FUNCTION Read_Bathy_level
220  !
221  !*****************************************************
222  !   function Write_Bathy_level(name,Grid)
223  !*****************************************************
224  !
225  INTEGER FUNCTION Write_Bathy_level(name,Grid)
226    !
227    USE io_netcdf
228    !     
229    CHARACTER(*) name
230    TYPE(Coordinates) :: Grid
231    INTEGER :: status,ncid
232    CHARACTER(len=1),DIMENSION(2) :: dimnames
233    !
234    status = nf90_create(name,NF90_WRITE,ncid)
235    status = nf90_close(ncid)         
236    !
237    dimnames = (/ 'x','y' /)
238    CALL Write_Ncdf_dim(dimnames(1),name,nxfin)
239    CALL Write_Ncdf_dim(dimnames(2),name,nyfin)
240    !     
241    CALL Write_Ncdf_var('nav_lon',dimnames,name,Grid%nav_lon    ,'float')
242    CALL Write_Ncdf_var('nav_lat',dimnames,name,Grid%nav_lat    ,'float')
243    CALL Write_Ncdf_var('mbathy' ,dimnames,name,Grid%bathy_level,'float')
244    !
245    CALL Copy_Ncdf_att('nav_lon',TRIM(parent_meshmask_file),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon))
246    CALL Copy_Ncdf_att('nav_lat',TRIM(parent_meshmask_file),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat))
247    CALL Copy_Ncdf_att('mbathy' ,TRIM(parent_meshmask_file),name)       
248    !
249    WRITE(*,*) ' '
250    WRITE(*,*) 'Writing bathymetry file: ',name
251    WRITE(*,*) ' '
252    !
253    Write_Bathy_level = 1
254    !
255  END FUNCTION Write_Bathy_level
256  !
257  !*****************************************************
258  !   function Read_Bathy_meter(name,CoarseGrid,ChildGrid)
259  !*****************************************************
260  !
261  INTEGER FUNCTION Read_Bathy_meter(name,CoarseGrid,ChildGrid,Pacifique)
262    !
263    USE io_netcdf
264    CHARACTER(*) name
265    INTEGER :: i,j,tabdim1,tabdim2
266    INTEGER, DIMENSION(1) :: i_min,i_max,j_min,j_max
267    REAL*8,POINTER,DIMENSION(:) :: topo_lon,topo_lat
268    INTEGER :: status,ncid,varid 
269    LOGICAL,OPTIONAL :: Pacifique
270    TYPE(Coordinates) :: CoarseGrid,ChildGrid
271    !
272    IF( Dims_Existence('lon',name) .AND. Dims_Existence('lat',name) ) THEN
273       WRITE(*,*) '****'
274       WRITE(*,*) ' etopo format for external high resolution database  '
275       WRITE(*,*) '****'
276       CALL Read_Ncdf_var('lon',name,topo_lon)
277       CALL Read_Ncdf_var('lat',name,topo_lat)
278    ELSE IF( Dims_Existence('x',name) .AND. Dims_Existence('y',name) ) THEN
279       WRITE(*,*) '****'
280       WRITE(*,*) ' OPA format for external high resolution database  '
281       WRITE(*,*) '****'
282       CALL Read_Ncdf_var('nav_lon',name,CoarseGrid%nav_lon)
283       CALL Read_Ncdf_var('nav_lat',name,CoarseGrid%nav_lat)
284       CALL Read_Ncdf_var(parent_batmet_name,name,CoarseGrid%Bathy_meter)
285       !           
286       IF ( PRESENT(Pacifique) ) THEN
287          IF(Pacifique) THEN
288             WHERE(CoarseGrid%nav_lon < 0.001) 
289                CoarseGrid%nav_lon = CoarseGrid%nav_lon + 360.
290             END WHERE
291          ENDIF
292       ENDIF
293       !     
294       Read_Bathy_meter = 1
295       RETURN     
296    ELSE
297       WRITE(*,*) '****'
298       WRITE(*,*) '*** ERROR Bad format for external high resolution database'
299       WRITE(*,*) '****'
300       STOP 
301    ENDIF
302    !
303    IF( MAXVAL(ChildGrid%glamt) > 180 ) THEN                 
304       !         
305       WHERE( topo_lon < 0 )
306          topo_lon = topo_lon + 360.
307       END WHERE
308       !         
309       i_min = MAXLOC(topo_lon,mask = topo_lon < MINVAL(ChildGrid%nav_lon))
310       i_max = MINLOC(topo_lon,mask = topo_lon > MAXVAL(ChildGrid%nav_lon))                   
311       j_min = MAXLOC(topo_lat,mask = topo_lat < MINVAL(ChildGrid%nav_lat))
312       j_max = MINLOC(topo_lat,mask = topo_lat > MAXVAL(ChildGrid%nav_lat))
313       !         
314       tabdim1 = ( SIZE(topo_lon) - i_min(1) + 1 ) + i_max(1)                   
315       !         
316       IF(j_min(1)-2 >= 1 .AND. j_max(1)+3 <= SIZE(topo_lat,1) ) THEN
317          j_min(1) = j_min(1)-2
318          j_max(1) = j_max(1)+3
319       ENDIF
320       tabdim2 = j_max(1) - j_min(1) + 1
321       !
322       ALLOCATE(CoarseGrid%nav_lon(tabdim1,tabdim2))
323       ALLOCATE(CoarseGrid%nav_lat(tabdim1,tabdim2))
324       ALLOCATE(CoarseGrid%Bathy_meter(tabdim1,tabdim2))         
325       !
326       DO i = 1,tabdim1
327          CoarseGrid%nav_lat(i,:) = topo_lat(j_min(1):j_max(1))
328       END DO
329       !
330       DO j = 1, tabdim2
331          !                     
332          CoarseGrid%nav_lon(1:SIZE(topo_lon)-i_min(1)+1      ,j) = topo_lon(i_min(1):SIZE(topo_lon))
333          CoarseGrid%nav_lon(2+SIZE(topo_lon)-i_min(1):tabdim1,j) = topo_lon(1:i_max(1))
334          !   
335       END DO
336       status = nf90_open(name,NF90_NOWRITE,ncid)
337       status = nf90_inq_varid(ncid,elevation_name,varid)
338       !         
339       status=nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter(1:SIZE(topo_lon)-i_min(1)+1,:), &
340            start=(/i_min(1),j_min(1)/),count=(/SIZE(topo_lon)-i_min(1),tabdim2/))
341
342       status=nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter(2+SIZE(topo_lon)-i_min(1):tabdim1,:), &
343            start=(/1,j_min(1)/),count=(/i_max(1),tabdim2/))               
344       !
345    ELSE
346       !
347       i_min = MAXLOC(topo_lon,mask = topo_lon < MINVAL(ChildGrid%nav_lon))
348       i_max = MINLOC(topo_lon,mask = topo_lon > MAXVAL(ChildGrid%nav_lon))
349       j_min = MAXLOC(topo_lat,mask = topo_lat < MINVAL(ChildGrid%nav_lat))
350       j_max = MINLOC(topo_lat,mask = topo_lat > MAXVAL(ChildGrid%nav_lat))
351       !     
352       IF(i_min(1)-2 >= 1 .AND. i_max(1)+3 <= SIZE(topo_lon,1) ) THEN
353          i_min(1) = i_min(1)-2
354          i_max(1) = i_max(1)+3
355       ENDIF
356       tabdim1 = i_max(1) - i_min(1) + 1
357       !
358       IF(j_min(1)-2 >= 1 .AND. j_max(1)+3 <= SIZE(topo_lat,1) ) THEN
359          j_min(1) = j_min(1)-2
360          j_max(1) = j_max(1)+3
361       ENDIF
362       tabdim2 = j_max(1) - j_min(1) + 1
363       !
364       WRITE(*,*) ' '
365       WRITE(*,*) 'Reading bathymetry file: ',name
366       WRITE(*,*) ' '
367       !
368       ALLOCATE(CoarseGrid%nav_lon(tabdim1,tabdim2))
369       ALLOCATE(CoarseGrid%nav_lat(tabdim1,tabdim2))
370       ALLOCATE(CoarseGrid%Bathy_meter(tabdim1,tabdim2))
371       !
372       DO j = 1,tabdim2
373          CoarseGrid%nav_lon(:,j) = topo_lon(i_min(1):i_max(1))
374       END DO
375       !     
376       DO i = 1,tabdim1
377          CoarseGrid%nav_lat(i,:) = topo_lat(j_min(1):j_max(1)) 
378       END DO
379       !
380       status = nf90_open(name,NF90_NOWRITE,ncid)
381       status = nf90_inq_varid(ncid,elevation_name,varid)
382       status = nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter, &
383          &                  start=(/i_min(1),j_min(1)/),count=(/tabdim1,tabdim2/))
384       !
385    ENDIF
386    !
387    status = nf90_close(ncid)     
388    !
389    WHERE(CoarseGrid%Bathy_meter.GE.0)
390       CoarseGrid%Bathy_meter = 0.0
391    END WHERE
392    !     
393    CoarseGrid%Bathy_meter(:,:) = -1.0 * CoarseGrid%Bathy_meter(:,:)
394    !     
395    Read_Bathy_meter = 1
396    RETURN
397    !     
398  END FUNCTION Read_Bathy_meter
399  !
400  !
401  !*****************************************************
402  !   function Read_Bathy_meter(name,CoarseGrid,ChildGrid)
403  !*****************************************************
404  !
405  INTEGER FUNCTION Read_Bathymeter(name,Grid)
406    !
407    !
408    USE io_netcdf
409    !     
410    CHARACTER(*) name
411    TYPE(Coordinates) :: Grid
412    !
413    CALL Read_Ncdf_var(parent_batmet_name,name,Grid%Bathy_meter)   
414    !
415    WRITE(*,*) ' '
416    WRITE(*,*) 'Reading bathymetry file: ',name
417    WRITE(*,*) ' '     
418    !
419    Read_Bathymeter = 1
420    !     
421  END FUNCTION Read_Bathymeter
422  !     
423  !*****************************************************
424  !   function Write_Bathy_meter(name,Grid)
425  !*****************************************************
426  !
427  INTEGER FUNCTION Write_Bathy_meter(name,Grid)
428    !
429    USE io_netcdf
430    !
431    CHARACTER(*) name
432    TYPE(Coordinates) :: Grid
433    INTEGER :: status,ncid
434    CHARACTER(len=1),DIMENSION(2) :: dimnames
435    INTEGER :: nx,ny
436    !
437    status = nf90_create(name,NF90_WRITE,ncid)
438    status = nf90_close(ncid)     
439    !
440    nx = SIZE(Grid%bathy_meter,1)
441    ny = SIZE(Grid%bathy_meter,2)
442    dimnames = (/ 'x','y' /)
443
444    CALL Write_Ncdf_dim(dimnames(1),name,nx)
445    CALL Write_Ncdf_dim(dimnames(2),name,ny)
446    !     
447    CALL Write_Ncdf_var('nav_lon'         ,dimnames,name,Grid%nav_lon    ,'float')
448    CALL Write_Ncdf_var('nav_lat'         ,dimnames,name,Grid%nav_lat    ,'float')
449    CALL Write_Ncdf_var(parent_batmet_name,dimnames,name,Grid%bathy_meter,'float')
450    !
451    CALL Copy_Ncdf_att('nav_lon'         ,TRIM(parent_bathy_meter),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon))
452    CALL Copy_Ncdf_att('nav_lat'         ,TRIM(parent_bathy_meter),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat)) 
453    CALL Copy_Ncdf_att(parent_batmet_name,TRIM(parent_bathy_meter),name)
454    !
455    WRITE(*,*) ' '
456    WRITE(*,*) 'Writing bathymetry file: ',name
457    WRITE(*,*) ' '
458    !
459    Write_Bathy_meter = 1
460    !     
461  END FUNCTION Write_Bathy_meter
462  !
463  !*****************************************************
464  !   function set_child_name(Parentname,Childname)
465  !*****************************************************
466  !
467  SUBROUTINE set_child_name(Parentname,Childname)
468    !
469    CHARACTER(*),INTENT(in) :: Parentname
470    CHARACTER(*),INTENT(out) :: Childname
471    CHARACTER(2) :: prefix
472    INTEGER :: pos
473    !   
474    pos  = INDEX(TRIM(Parentname),'/',back=.TRUE.)
475    !
476    prefix=Parentname(pos+1:pos+2)
477    IF (prefix == '1_') THEN
478       Childname = '2_'//Parentname(pos+3:LEN(Parentname)) 
479    ELSEIF (prefix == '2_') THEN
480       Childname = '3_'//Parentname(pos+3:LEN(Parentname)) 
481    ELSEIF (prefix == '3_') THEN
482       Childname = '4_'//Parentname(pos+3:LEN(Parentname)) 
483    ELSEIF (prefix == '4_') THEN
484       Childname = '5_'//Parentname(pos+3:LEN(Parentname)) 
485    ELSE
486       Childname = '1_'//Parentname(pos+1:LEN(Parentname)) 
487    ENDIF
488    !   
489  END SUBROUTINE set_child_name
490  !
491  !*****************************************************
492  !   function set_child_name(Parentname,Childname)
493  !*****************************************************
494  !
495  !*****************************************************
496  !   subroutine get_interptype(varname,interp_type,conservation)
497  !*****************************************************
498  !
499  SUBROUTINE get_interptype( varname,interp_type,conservation )
500    !
501    LOGICAL,OPTIONAL :: conservation
502    CHARACTER(*) :: interp_type,varname
503    INTEGER :: pos,pos1,pos2,pos3,i,k 
504    LOGICAL :: find       
505    i=1
506    DO WHILE ( TRIM(VAR_INTERP(i)) .NE. 'NULL' )     
507       pos = INDEX( TRIM(VAR_INTERP(i)) , TRIM(varname) )
508       IF ( pos .NE. 0 ) THEN     
509          pos1 = INDEX( TRIM(VAR_INTERP(i)) , 'bicubic' )
510          pos2 = INDEX( TRIM(VAR_INTERP(i)) , 'bilinear' )
511          pos3 = INDEX( TRIM(VAR_INTERP(i)) , 'conservative' )
512          ! initialize interp_type
513          IF( pos1 .NE. 0 ) interp_type = 'bicubic'
514          IF( pos2 .NE. 0 ) interp_type = 'bilinear'
515          IF( pos1 .EQ. 0 .AND. pos2 .EQ. 0) interp_type = 'bicubic'
516          ! initialize conservation                                           
517          IF( pos3 .NE. 0 .AND. PRESENT(conservation) ) THEN
518             conservation = .TRUE.
519             RETURN
520          ELSE
521             conservation = .FALSE.
522          ENDIF
523          find = .FALSE.
524          IF( PRESENT(conservation) ) THEN
525             k=0
526             conservation = .FALSE.         
527             DO WHILE( k < SIZE(flxtab) .AND. .NOT.find ) 
528                k = k+1
529                IF( TRIM(varname) .EQ. TRIM(flxtab(k)) ) THEN       
530                   conservation = .TRUE.
531                   find = .TRUE.
532                ENDIF
533             END DO
534          ENDIF
535          RETURN
536       ENDIF
537       i = i+1
538    END DO
539    !default values interp_type = bicubic // conservation = false     
540    interp_type = 'bicubic'     
541    IF( PRESENT(conservation) ) conservation = .FALSE.
542
543    RETURN
544    !   
545  END SUBROUTINE get_interptype
546  !     
547  !*****************************************************             
548  !   end subroutine get_interptype
549  !*****************************************************
550  !           
551  !*****************************************************
552  !   subroutine Init_mask(name,Grid)
553  !*****************************************************
554  !
555  SUBROUTINE Init_mask(name,Grid,jpiglo,jpjglo)
556    !
557    USE io_netcdf
558    !     
559    CHARACTER(*) name
560    INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo
561    TYPE(Coordinates) :: Grid
562    REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL()
563    !
564    IF(jpiglo == 1 .AND. jpjglo == 1) THEN
565       CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level)
566    ELSE
567       CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) )
568    ENDIF
569
570
571    !
572    WRITE(*,*) 'Init masks in T,U,V,F points'   
573    !
574    nx = SIZE(Grid%Bathy_level,1)
575    ny = SIZE(Grid%Bathy_level,2)
576    !
577    !     
578    ALLOCATE(Grid%tmask(nx,ny,N), &
579         Grid%umask(nx,ny,N), &
580         Grid%vmask(nx,ny,N), &
581         Grid%fmask(nx,ny,N))
582    !
583    DO k = 1,N
584       !     
585       WHERE(Grid%Bathy_level(:,:) <= k-1 )   
586          Grid%tmask(:,:,k) = 0
587       ELSEWHERE
588          Grid%tmask(:,:,k) = 1
589       END WHERE
590       !
591    END DO
592    !
593    Grid%umask(1:nx-1,:,:) = Grid%tmask(1:nx-1,:,:)*Grid%tmask(2:nx,:,:)
594    Grid%vmask(:,1:ny-1,:) = Grid%tmask(:,1:ny-1,:)*Grid%tmask(:,2:ny,:)
595    !
596    Grid%umask(nx,:,:) = Grid%tmask(nx,:,:)
597    Grid%vmask(:,ny,:) = Grid%tmask(:,ny,:)
598    !     
599    Grid%fmask(1:nx-1,1:ny-1,:) = Grid%tmask(1:nx-1,1:ny-1,:)*Grid%tmask(2:nx,1:ny-1,:)* &
600         Grid%tmask(1:nx-1,2:ny,:)*Grid%tmask(2:nx,2:ny,:) 
601    !
602    Grid%fmask(nx,:,:) = Grid%tmask(nx,:,:)
603    Grid%fmask(:,ny,:) = Grid%tmask(:,ny,:)
604    !
605    ALLOCATE(zwf(nx,ny))
606    !     
607    DO k = 1,N
608       !
609       zwf(:,:) = Grid%fmask(:,:,k)
610       !         
611       DO j = 2, ny-1
612          DO i = 2,nx-1           
613             IF( Grid%fmask(i,j,k) == 0. ) THEN               
614                Grid%fmask(i,j,k) = shlat * MIN(1.,MAX( zwf(i+1,j),zwf(i,j+1),zwf(i-1,j),zwf(i,j-1)))
615             END IF
616          END DO
617       END DO
618       !
619       DO j = 2, ny-1
620          IF( Grid%fmask(1,j,k) == 0. ) THEN
621             Grid%fmask(1,j,k) = shlat * MIN(1.,MAX(zwf(2,j),zwf(1,j+1),zwf(1,j-1)))
622          ENDIF
623
624          IF( Grid%fmask(nx,j,k) == 0. ) THEN
625             Grid%fmask(nx,j,k) = shlat * MIN(1.,MAX(zwf(nx,j+1),zwf(nx-1,j),zwf(nx,j-1)))
626          ENDIF
627       END DO
628       !         
629       DO i = 2, nx-1         
630          IF( Grid%fmask(i,1,k) == 0. ) THEN
631             Grid%fmask(i, 1 ,k) = shlat*MIN( 1.,MAX(zwf(i+1,1),zwf(i,2),zwf(i-1,1)))
632          ENDIF
633          !           
634          IF( Grid%fmask(i,ny,k) == 0. ) THEN
635             Grid%fmask(i,ny,k) = shlat * MIN(1.,MAX(zwf(i+1,ny),zwf(i-1,ny),zwf(i,ny-1)))
636          ENDIF
637       END DO
638       !!
639    END DO
640    !!     
641  END SUBROUTINE Init_mask
642  !
643  !*****************************************************
644  !   end subroutine Init_mask
645  !*****************************************************
646  !
647  !*****************************************************
648  !   subroutine Init_Tmask(name,Grid)
649  !*****************************************************
650  !
651  SUBROUTINE Init_Tmask(name,Grid,jpiglo,jpjglo)
652    !
653    USE io_netcdf
654    !     
655    CHARACTER(*) name
656    INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo
657    TYPE(Coordinates) :: Grid
658    REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL()
659    !
660    IF(jpiglo == 1 .AND. jpjglo == 1) THEN
661       CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level)
662    ELSE
663       CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) )
664    ENDIF
665    !
666    nx = SIZE(Grid%Bathy_level,1)
667    ny = SIZE(Grid%Bathy_level,2) 
668    !
669    WRITE(*,*) 'Init masks in T points'   
670    !     
671    ALLOCATE(Grid%tmask(nx,ny,N))
672    !
673    DO k = 1,N
674       !     
675       WHERE(Grid%Bathy_level(:,:) <= k-1 )   
676          Grid%tmask(:,:,k) = 0.
677       ELSEWHERE
678          Grid%tmask(:,:,k) = 1.
679       END WHERE
680       !
681    END DO
682    !     
683  END SUBROUTINE Init_Tmask
684  !
685  !*****************************************************
686  !   subroutine get_mask(name,Grid)
687  !*****************************************************
688  !
689  SUBROUTINE get_mask(level,posvar,mask,filename)
690    !
691    USE io_netcdf
692    !     
693    CHARACTER(*) filename
694    CHARACTER(*) posvar
695    INTEGER :: level, nx, ny
696    LOGICAL,DIMENSION(:,:),POINTER :: mask
697    INTEGER,DIMENSION(:,:),POINTER :: maskT,maskU,maskV
698    !     
699    TYPE(Coordinates) :: Grid
700    !
701    CALL Read_Ncdf_var('Bathy_level',filename,Grid%Bathy_level)
702    !
703    nx = SIZE(Grid%Bathy_level,1)
704    ny = SIZE(Grid%Bathy_level,2)
705    ALLOCATE(maskT(nx,ny),mask(nx,ny))
706    mask = .TRUE.
707    !
708    WHERE(Grid%Bathy_level(:,:) <= level-1 )   
709       maskT(:,:) = 0
710    ELSEWHERE
711       maskT(:,:) = 1
712    END WHERE
713    !
714    SELECT CASE(posvar)
715       !
716    CASE('T')
717       !
718       WHERE(maskT > 0)
719          mask = .TRUE.
720       ELSEWHERE
721          mask = .FALSE.
722       END WHERE
723       DEALLOCATE(maskT)
724       !
725    CASE('U')
726       !
727       ALLOCATE(maskU(nx,ny))
728       maskU(1:nx-1,:) = maskT(1:nx-1,:)*maskT(2:nx,:)
729       maskU(nx,:) = maskT(nx,:)
730       WHERE(maskU > 0)
731          mask = .TRUE.
732       ELSEWHERE
733          mask = .FALSE.
734       END WHERE
735       DEALLOCATE(maskU,maskT)
736       !
737    CASE('V')
738       !
739       ALLOCATE(maskV(nx,ny))
740       maskV(:,1:ny-1) = maskT(:,1:ny-1)*maskT(:,2:ny)
741       maskV(:,ny) = maskT(:,ny)
742       WHERE(maskT > 0)
743          mask = .TRUE.
744       ELSEWHERE
745          mask = .FALSE.
746       END WHERE
747       DEALLOCATE(maskV,maskT)
748       !
749    END SELECT
750    !     
751  END SUBROUTINE get_mask
752  !
753  !*****************************************************
754  !   end subroutine get_mask
755  !*****************************************************
756  !       
757  !
758  !*****************************************************
759  !   subroutine read_dimg_var(unit,irec,field)
760  !*****************************************************
761  !
762  SUBROUTINE read_dimg_var(unit,irec,field,jpk)
763    !     
764    INTEGER :: unit,irec,jpk
765    REAL*8,DIMENSION(:,:,:,:),POINTER :: field
766    INTEGER :: k
767    !     
768    DO k = 1,jpk
769       READ(unit,REC=irec) field(:,:,k,1)
770       irec = irec + 1
771    ENDDO
772    !
773  END SUBROUTINE read_dimg_var
774  !
775  !
776  !*****************************************************
777  !   subroutine read_dimg_var(unit,irec,field)
778  !*****************************************************
779  !
780  SUBROUTINE write_dimg_var(unit,irec,field,jpk)
781    !     
782    INTEGER :: unit,irec,jpk
783    REAL*8,DIMENSION(:,:,:,:),POINTER :: field
784    INTEGER :: k
785    !     
786    DO k = 1,jpk
787       WRITE(unit,REC=irec) field(:,:,k,1)
788       irec = irec + 1
789    ENDDO
790    !
791  END SUBROUTINE write_dimg_var
792!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!             
793!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
794END MODULE agrif_readwrite
Note: See TracBrowser for help on using the repository browser.