New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
agrif_readwrite.f90 in utils/tools/NESTING/src – NEMO

source: utils/tools/NESTING/src/agrif_readwrite.f90 @ 9628

Last change on this file since 9628 was 9628, checked in by jchanut, 7 years ago

Remove bathymetry connection if open boundary is closed, #2089
  • Property svn:keywords set to Id
File size: 26.9 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_NOCLOBBER,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    REAL*8 :: zdel
272    zdel = 0.5 ! Offset in degrees to extend extraction of bathymetry data
273    !
274    IF( Dims_Existence('lon',name) .AND. Dims_Existence('lat',name) ) THEN
275       WRITE(*,*) '****'
276       WRITE(*,*) ' etopo format for external high resolution database  '
277       WRITE(*,*) '****'
278       CALL Read_Ncdf_var('lon',name,topo_lon)
279       CALL Read_Ncdf_var('lat',name,topo_lat)
280    ELSE IF( Dims_Existence('x',name) .AND. Dims_Existence('y',name) ) THEN
281       WRITE(*,*) '****'
282       WRITE(*,*) ' OPA format for external high resolution database  '
283       WRITE(*,*) '****'
284       CALL Read_Ncdf_var('nav_lon',name,CoarseGrid%nav_lon)
285       CALL Read_Ncdf_var('nav_lat',name,CoarseGrid%nav_lat)
286       CALL Read_Ncdf_var(parent_batmet_name,name,CoarseGrid%Bathy_meter)
287       !           
288       IF ( PRESENT(Pacifique) ) THEN
289          IF(Pacifique) THEN
290             WHERE(CoarseGrid%nav_lon < 0.001) 
291                CoarseGrid%nav_lon = CoarseGrid%nav_lon + 360.
292             END WHERE
293          ENDIF
294       ENDIF
295       !     
296       Read_Bathy_meter = 1
297       RETURN     
298    ELSE
299       WRITE(*,*) '****'
300       WRITE(*,*) '*** ERROR Bad format for external high resolution database'
301       WRITE(*,*) '****'
302       STOP 
303    ENDIF
304    !
305    IF( MAXVAL(ChildGrid%glamt) > 180 ) THEN                 
306       !         
307       WHERE( topo_lon < 0 )
308          topo_lon = topo_lon + 360.
309       END WHERE
310       !         
311       i_min = MAXLOC(topo_lon,mask = topo_lon < MINVAL(ChildGrid%nav_lon)-zdel)
312       i_max = MINLOC(topo_lon,mask = topo_lon > MAXVAL(ChildGrid%nav_lon)+zdel)                   
313       j_min = MAXLOC(topo_lat,mask = topo_lat < MINVAL(ChildGrid%nav_lat)-zdel)
314       j_max = MINLOC(topo_lat,mask = topo_lat > MAXVAL(ChildGrid%nav_lat)+zdel)
315       !         
316       tabdim1 = ( SIZE(topo_lon) - i_min(1) + 1 ) + i_max(1)                   
317       !         
318       IF(j_min(1)-2 >= 1 .AND. j_max(1)+3 <= SIZE(topo_lat,1) ) THEN
319          j_min(1) = j_min(1)-2
320          j_max(1) = j_max(1)+3
321       ENDIF
322       tabdim2 = j_max(1) - j_min(1) + 1
323       !
324       ALLOCATE(CoarseGrid%nav_lon(tabdim1,tabdim2))
325       ALLOCATE(CoarseGrid%nav_lat(tabdim1,tabdim2))
326       ALLOCATE(CoarseGrid%Bathy_meter(tabdim1,tabdim2))         
327       !
328       DO i = 1,tabdim1
329          CoarseGrid%nav_lat(i,:) = topo_lat(j_min(1):j_max(1))
330       END DO
331       !
332       DO j = 1, tabdim2
333          !                     
334          CoarseGrid%nav_lon(1:SIZE(topo_lon)-i_min(1)+1      ,j) = topo_lon(i_min(1):SIZE(topo_lon))
335          CoarseGrid%nav_lon(2+SIZE(topo_lon)-i_min(1):tabdim1,j) = topo_lon(1:i_max(1))
336          !   
337       END DO
338       status = nf90_open(name,NF90_NOWRITE,ncid)
339       status = nf90_inq_varid(ncid,elevation_name,varid)
340       !
341       IF (status/=nf90_noerr) THEN
342          WRITE(*,*)"Can't find variable: ", elevation_name
343          STOP
344       ENDIF
345       !         
346       status=nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter(1:SIZE(topo_lon)-i_min(1)+1,:), &
347            start=(/i_min(1),j_min(1)/),count=(/SIZE(topo_lon)-i_min(1),tabdim2/))
348
349       status=nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter(2+SIZE(topo_lon)-i_min(1):tabdim1,:), &
350            start=(/1,j_min(1)/),count=(/i_max(1),tabdim2/))               
351       !
352    ELSE
353       !
354       i_min = MAXLOC(topo_lon,mask = topo_lon < MINVAL(ChildGrid%nav_lon)-zdel)
355       i_max = MINLOC(topo_lon,mask = topo_lon > MAXVAL(ChildGrid%nav_lon)+zdel)
356       j_min = MAXLOC(topo_lat,mask = topo_lat < MINVAL(ChildGrid%nav_lat)-zdel)
357       j_max = MINLOC(topo_lat,mask = topo_lat > MAXVAL(ChildGrid%nav_lat)+zdel)
358       !     
359       IF(i_min(1)-2 >= 1 .AND. i_max(1)+3 <= SIZE(topo_lon,1) ) THEN
360          i_min(1) = i_min(1)-2
361          i_max(1) = i_max(1)+3
362       ENDIF
363       tabdim1 = i_max(1) - i_min(1) + 1
364       !
365       IF(j_min(1)-2 >= 1 .AND. j_max(1)+3 <= SIZE(topo_lat,1) ) THEN
366          j_min(1) = j_min(1)-2
367          j_max(1) = j_max(1)+3
368       ENDIF
369       tabdim2 = j_max(1) - j_min(1) + 1
370       !
371       WRITE(*,*) ' '
372       WRITE(*,*) 'Reading bathymetry file: ',name
373       WRITE(*,*) ' '
374       !
375       ALLOCATE(CoarseGrid%nav_lon(tabdim1,tabdim2))
376       ALLOCATE(CoarseGrid%nav_lat(tabdim1,tabdim2))
377       ALLOCATE(CoarseGrid%Bathy_meter(tabdim1,tabdim2))
378       !
379       DO j = 1,tabdim2
380          CoarseGrid%nav_lon(:,j) = topo_lon(i_min(1):i_max(1))
381       END DO
382       !     
383       DO i = 1,tabdim1
384          CoarseGrid%nav_lat(i,:) = topo_lat(j_min(1):j_max(1)) 
385       END DO
386       !
387       status = nf90_open(name,NF90_NOWRITE,ncid)
388       status = nf90_inq_varid(ncid,elevation_name,varid)
389       IF (status/=nf90_noerr) THEN
390          WRITE(*,*)"Can't find variable: ", elevation_name
391          STOP
392       ENDIF
393
394       status = nf90_get_var(ncid,varid,CoarseGrid%Bathy_meter, &
395          &                  start=(/i_min(1),j_min(1)/),count=(/tabdim1,tabdim2/))
396       !
397    ENDIF
398    !
399    status = nf90_close(ncid)     
400    !
401    WHERE(CoarseGrid%Bathy_meter.GE.0)
402       CoarseGrid%Bathy_meter = 0.0
403    END WHERE
404    !     
405    CoarseGrid%Bathy_meter(:,:) = -1.0 * CoarseGrid%Bathy_meter(:,:)
406    !     
407    Read_Bathy_meter = 1
408    RETURN
409    !     
410  END FUNCTION Read_Bathy_meter
411  !
412  !
413  !*****************************************************
414  !   function Read_Bathy_meter(name,CoarseGrid,ChildGrid)
415  !*****************************************************
416  !
417  INTEGER FUNCTION Read_Bathymeter(name,Grid)
418    !
419    !
420    USE io_netcdf
421    !     
422    CHARACTER(*) name
423    TYPE(Coordinates) :: Grid
424    !
425    CALL Read_Ncdf_var(parent_batmet_name,name,Grid%Bathy_meter)   
426    !
427    WRITE(*,*) ' '
428    WRITE(*,*) 'Reading bathymetry file: ',name
429    WRITE(*,*) ' '     
430    !
431    Read_Bathymeter = 1
432    !     
433  END FUNCTION Read_Bathymeter
434  !     
435  !*****************************************************
436  !   function Write_Bathy_meter(name,Grid)
437  !*****************************************************
438  !
439  INTEGER FUNCTION Write_Bathy_meter(name,Grid)
440    !
441    USE io_netcdf
442    !
443    CHARACTER(*) name
444    TYPE(Coordinates) :: Grid
445    INTEGER :: status,ncid
446    CHARACTER(len=1),DIMENSION(2) :: dimnames
447    INTEGER :: nx,ny
448    !
449    status = nf90_create(name,NF90_WRITE,ncid)
450    status = nf90_close(ncid)     
451    !
452    nx = SIZE(Grid%bathy_meter,1)
453    ny = SIZE(Grid%bathy_meter,2)
454    dimnames = (/ 'x','y' /)
455
456    CALL Write_Ncdf_dim(dimnames(1),name,nx)
457    CALL Write_Ncdf_dim(dimnames(2),name,ny)
458    !     
459    CALL Write_Ncdf_var('nav_lon'         ,dimnames,name,Grid%nav_lon    ,'float')
460    CALL Write_Ncdf_var('nav_lat'         ,dimnames,name,Grid%nav_lat    ,'float')
461    CALL Write_Ncdf_var(parent_batmet_name,dimnames,name,Grid%bathy_meter,'float')
462    CALL Write_Ncdf_var('weight'          ,dimnames,name,Grid%wgt        ,'float')
463    !
464    CALL Copy_Ncdf_att('nav_lon'         ,TRIM(parent_bathy_meter),name,MINVAL(Grid%nav_lon),MAXVAL(Grid%nav_lon))
465    CALL Copy_Ncdf_att('nav_lat'         ,TRIM(parent_bathy_meter),name,MINVAL(Grid%nav_lat),MAXVAL(Grid%nav_lat)) 
466    CALL Copy_Ncdf_att(parent_batmet_name,TRIM(parent_bathy_meter),name)
467    !
468    WRITE(*,*) ' '
469    WRITE(*,*) 'Writing bathymetry file: ',name
470    WRITE(*,*) ' '
471    !
472    Write_Bathy_meter = 1
473    !     
474  END FUNCTION Write_Bathy_meter
475  !
476  !*****************************************************
477  !   function set_child_name(Parentname,Childname)
478  !*****************************************************
479  !
480  SUBROUTINE set_child_name(Parentname,Childname)
481    !
482    CHARACTER(*),INTENT(in) :: Parentname
483    CHARACTER(*),INTENT(out) :: Childname
484    CHARACTER(2) :: prefix
485    INTEGER :: pos
486    !   
487    pos  = INDEX(TRIM(Parentname),'/',back=.TRUE.)
488    !
489    prefix=Parentname(pos+1:pos+2)
490    IF (prefix == '1_') THEN
491       Childname = '2_'//Parentname(pos+3:LEN(Parentname)) 
492    ELSEIF (prefix == '2_') THEN
493       Childname = '3_'//Parentname(pos+3:LEN(Parentname)) 
494    ELSEIF (prefix == '3_') THEN
495       Childname = '4_'//Parentname(pos+3:LEN(Parentname)) 
496    ELSEIF (prefix == '4_') THEN
497       Childname = '5_'//Parentname(pos+3:LEN(Parentname)) 
498    ELSE
499       Childname = '1_'//Parentname(pos+1:LEN(Parentname)) 
500    ENDIF
501    !   
502  END SUBROUTINE set_child_name
503  !
504  !*****************************************************
505  !   function set_child_name(Parentname,Childname)
506  !*****************************************************
507  !
508  !*****************************************************
509  !   subroutine get_interptype(varname,interp_type,conservation)
510  !*****************************************************
511  !
512  SUBROUTINE get_interptype( varname,interp_type,conservation )
513    !
514    LOGICAL,OPTIONAL :: conservation
515    CHARACTER(*) :: interp_type,varname
516    INTEGER :: pos,pos1,pos2,pos3,i,k 
517    LOGICAL :: find       
518    i=1
519    DO WHILE ( TRIM(VAR_INTERP(i)) .NE. 'NULL' )     
520       pos = INDEX( TRIM(VAR_INTERP(i)) , TRIM(varname) )
521       IF ( pos .NE. 0 ) THEN     
522          pos1 = INDEX( TRIM(VAR_INTERP(i)) , 'bicubic' )
523          pos2 = INDEX( TRIM(VAR_INTERP(i)) , 'bilinear' )
524          pos3 = INDEX( TRIM(VAR_INTERP(i)) , 'conservative' )
525          ! initialize interp_type
526          IF( pos1 .NE. 0 ) interp_type = 'bicubic'
527          IF( pos2 .NE. 0 ) interp_type = 'bilinear'
528          IF( pos1 .EQ. 0 .AND. pos2 .EQ. 0) interp_type = 'bicubic'
529          ! initialize conservation                                           
530          IF( pos3 .NE. 0 .AND. PRESENT(conservation) ) THEN
531             conservation = .TRUE.
532             RETURN
533          ELSE
534             conservation = .FALSE.
535          ENDIF
536          find = .FALSE.
537          IF( PRESENT(conservation) ) THEN
538             k=0
539             conservation = .FALSE.         
540             DO WHILE( k < SIZE(flxtab) .AND. .NOT.find ) 
541                k = k+1
542                IF( TRIM(varname) .EQ. TRIM(flxtab(k)) ) THEN       
543                   conservation = .TRUE.
544                   find = .TRUE.
545                ENDIF
546             END DO
547          ENDIF
548          RETURN
549       ENDIF
550       i = i+1
551    END DO
552    !default values interp_type = bicubic // conservation = false     
553    interp_type = 'bicubic'     
554    IF( PRESENT(conservation) ) conservation = .FALSE.
555
556    RETURN
557    !   
558  END SUBROUTINE get_interptype
559  !     
560  !*****************************************************             
561  !   end subroutine get_interptype
562  !*****************************************************
563  !           
564  !*****************************************************
565  !   subroutine Init_mask(name,Grid)
566  !*****************************************************
567  !
568  SUBROUTINE Init_mask(name,Grid,jpiglo,jpjglo)
569    !
570    USE io_netcdf
571    !     
572    CHARACTER(*) name
573    INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo
574    TYPE(Coordinates) :: Grid
575    REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL()
576    !
577    IF(jpiglo == 1 .AND. jpjglo == 1) THEN
578       CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level)
579    ELSE
580       CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) )
581    ENDIF
582
583
584    !
585    WRITE(*,*) 'Init masks in T,U,V,F points'   
586    !
587    nx = SIZE(Grid%Bathy_level,1)
588    ny = SIZE(Grid%Bathy_level,2)
589    !
590    !     
591    ALLOCATE(Grid%tmask(nx,ny,N), &
592         Grid%umask(nx,ny,N), &
593         Grid%vmask(nx,ny,N), &
594         Grid%fmask(nx,ny,N))
595    !
596    DO k = 1,N
597       !     
598       WHERE(Grid%Bathy_level(:,:) <= k-1 )   
599          Grid%tmask(:,:,k) = 0
600       ELSEWHERE
601          Grid%tmask(:,:,k) = 1
602       END WHERE
603       !
604    END DO
605    !
606    Grid%umask(1:nx-1,:,:) = Grid%tmask(1:nx-1,:,:)*Grid%tmask(2:nx,:,:)
607    Grid%vmask(:,1:ny-1,:) = Grid%tmask(:,1:ny-1,:)*Grid%tmask(:,2:ny,:)
608    !
609    Grid%umask(nx,:,:) = Grid%tmask(nx,:,:)
610    Grid%vmask(:,ny,:) = Grid%tmask(:,ny,:)
611    !     
612    Grid%fmask(1:nx-1,1:ny-1,:) = Grid%tmask(1:nx-1,1:ny-1,:)*Grid%tmask(2:nx,1:ny-1,:)* &
613         Grid%tmask(1:nx-1,2:ny,:)*Grid%tmask(2:nx,2:ny,:) 
614    !
615    Grid%fmask(nx,:,:) = Grid%tmask(nx,:,:)
616    Grid%fmask(:,ny,:) = Grid%tmask(:,ny,:)
617    !
618    ALLOCATE(zwf(nx,ny))
619    !     
620    DO k = 1,N
621       !
622       zwf(:,:) = Grid%fmask(:,:,k)
623       !         
624       DO j = 2, ny-1
625          DO i = 2,nx-1           
626             IF( Grid%fmask(i,j,k) == 0. ) THEN               
627                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)))
628             END IF
629          END DO
630       END DO
631       !
632       DO j = 2, ny-1
633          IF( Grid%fmask(1,j,k) == 0. ) THEN
634             Grid%fmask(1,j,k) = shlat * MIN(1.,MAX(zwf(2,j),zwf(1,j+1),zwf(1,j-1)))
635          ENDIF
636
637          IF( Grid%fmask(nx,j,k) == 0. ) THEN
638             Grid%fmask(nx,j,k) = shlat * MIN(1.,MAX(zwf(nx,j+1),zwf(nx-1,j),zwf(nx,j-1)))
639          ENDIF
640       END DO
641       !         
642       DO i = 2, nx-1         
643          IF( Grid%fmask(i,1,k) == 0. ) THEN
644             Grid%fmask(i, 1 ,k) = shlat*MIN( 1.,MAX(zwf(i+1,1),zwf(i,2),zwf(i-1,1)))
645          ENDIF
646          !           
647          IF( Grid%fmask(i,ny,k) == 0. ) THEN
648             Grid%fmask(i,ny,k) = shlat * MIN(1.,MAX(zwf(i+1,ny),zwf(i-1,ny),zwf(i,ny-1)))
649          ENDIF
650       END DO
651       !!
652    END DO
653    !!     
654  END SUBROUTINE Init_mask
655  !
656  !*****************************************************
657  !   end subroutine Init_mask
658  !*****************************************************
659  !
660  !*****************************************************
661  !   subroutine Init_Tmask(name,Grid)
662  !*****************************************************
663  !
664  SUBROUTINE Init_Tmask(name,Grid,jpiglo,jpjglo)
665    !
666    USE io_netcdf
667    !     
668    CHARACTER(*) name
669    INTEGER :: nx,ny,k,i,j,jpiglo,jpjglo
670    TYPE(Coordinates) :: Grid
671    REAL*8, POINTER, DIMENSION(:,:) ::zwf => NULL()
672    !
673    IF(jpiglo == 1 .AND. jpjglo == 1) THEN
674       CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level)
675    ELSE
676       CALL Read_Ncdf_var('Bathy_level',name,Grid%Bathy_level,(/jpizoom,jpjzoom/),(/jpiglo,jpjglo/) )
677    ENDIF
678    !
679    nx = SIZE(Grid%Bathy_level,1)
680    ny = SIZE(Grid%Bathy_level,2) 
681    !
682    WRITE(*,*) 'Init masks in T points'   
683    !     
684    ALLOCATE(Grid%tmask(nx,ny,N))
685    !
686    DO k = 1,N
687       !     
688       WHERE(Grid%Bathy_level(:,:) <= k-1 )   
689          Grid%tmask(:,:,k) = 0.
690       ELSEWHERE
691          Grid%tmask(:,:,k) = 1.
692       END WHERE
693       !
694    END DO
695    !     
696  END SUBROUTINE Init_Tmask
697  !
698  !*****************************************************
699  !   subroutine get_mask(name,Grid)
700  !*****************************************************
701  !
702  SUBROUTINE get_mask(level,posvar,mask,filename)
703    !
704    USE io_netcdf
705    !     
706    CHARACTER(*) filename
707    CHARACTER(*) posvar
708    INTEGER :: level, nx, ny
709    LOGICAL,DIMENSION(:,:),POINTER :: mask
710    INTEGER,DIMENSION(:,:),POINTER :: maskT,maskU,maskV
711    !     
712    TYPE(Coordinates) :: Grid
713    !
714    CALL Read_Ncdf_var('Bathy_level',filename,Grid%Bathy_level)
715    !
716    nx = SIZE(Grid%Bathy_level,1)
717    ny = SIZE(Grid%Bathy_level,2)
718    ALLOCATE(maskT(nx,ny),mask(nx,ny))
719    mask = .TRUE.
720    !
721    WHERE(Grid%Bathy_level(:,:) <= level-1 )   
722       maskT(:,:) = 0
723    ELSEWHERE
724       maskT(:,:) = 1
725    END WHERE
726    !
727    SELECT CASE(posvar)
728       !
729    CASE('T')
730       !
731       WHERE(maskT > 0)
732          mask = .TRUE.
733       ELSEWHERE
734          mask = .FALSE.
735       END WHERE
736       DEALLOCATE(maskT)
737       !
738    CASE('U')
739       !
740       ALLOCATE(maskU(nx,ny))
741       maskU(1:nx-1,:) = maskT(1:nx-1,:)*maskT(2:nx,:)
742       maskU(nx,:) = maskT(nx,:)
743       WHERE(maskU > 0)
744          mask = .TRUE.
745       ELSEWHERE
746          mask = .FALSE.
747       END WHERE
748       DEALLOCATE(maskU,maskT)
749       !
750    CASE('V')
751       !
752       ALLOCATE(maskV(nx,ny))
753       maskV(:,1:ny-1) = maskT(:,1:ny-1)*maskT(:,2:ny)
754       maskV(:,ny) = maskT(:,ny)
755       WHERE(maskT > 0)
756          mask = .TRUE.
757       ELSEWHERE
758          mask = .FALSE.
759       END WHERE
760       DEALLOCATE(maskV,maskT)
761       !
762    END SELECT
763    !     
764  END SUBROUTINE get_mask
765  !
766  !*****************************************************
767  !   end subroutine get_mask
768  !*****************************************************
769  !       
770  !
771  !*****************************************************
772  !   subroutine read_dimg_var(unit,irec,field)
773  !*****************************************************
774  !
775  SUBROUTINE read_dimg_var(unit,irec,field,jpk)
776    !     
777    INTEGER :: unit,irec,jpk
778    REAL*8,DIMENSION(:,:,:,:),POINTER :: field
779    INTEGER :: k
780    !     
781    DO k = 1,jpk
782       READ(unit,REC=irec) field(:,:,k,1)
783       irec = irec + 1
784    ENDDO
785    !
786  END SUBROUTINE read_dimg_var
787  !
788  !
789  !*****************************************************
790  !   subroutine read_dimg_var(unit,irec,field)
791  !*****************************************************
792  !
793  SUBROUTINE write_dimg_var(unit,irec,field,jpk)
794    !     
795    INTEGER :: unit,irec,jpk
796    REAL*8,DIMENSION(:,:,:,:),POINTER :: field
797    INTEGER :: k
798    !     
799    DO k = 1,jpk
800       WRITE(unit,REC=irec) field(:,:,k,1)
801       irec = irec + 1
802    ENDDO
803    !
804  END SUBROUTINE write_dimg_var
805!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!             
806!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
807END MODULE agrif_readwrite
Note: See TracBrowser for help on using the repository browser.