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_create_restart.f90 in utils/tools/NESTING/src – NEMO

source: utils/tools/NESTING/src/agrif_create_restart.f90 @ 10393

Last change on this file since 10393 was 10393, checked in by clem, 5 years ago

ocean restart from nesting tool should now work. BGC restart must still be debugged but it is out of my expertise
  • Property svn:keywords set to Id
File size: 21.6 KB
Line 
1!
2!************************************************************************
3! Fortran 95 OPA Nesting tools                  *
4!                          *
5!     Copyright (C) 2005 Florian Lemari�(Florian.Lemarie@imag.fr) *
6!                          *
7!************************************************************************
8!
9PROGRAM create_rstrt
10  !
11  USE NETCDF
12  USE bilinear_interp
13  USE bicubic_interp
14  USE agrif_readwrite
15  USE io_netcdf 
16  USE agrif_extrapolation
17  USE agrif_interpolation
18  USE agrif_partial_steps       
19  USE agrif_connect_topo
20  !
21  IMPLICIT NONE
22  !
23  !************************************************************************
24  !                           *
25  ! PROGRAM  CREATE_RSTRT                    *
26  !                           *
27  ! program to interpolate parent grid restart file to child grid    *
28  !                           *
29  !                           *
30  !Interpolation is carried out using bilinear interpolation      *
31  !routine from SCRIP package                *     
32  !                           *
33  !http://climate.lanl.gov/Software/SCRIP/            *
34  !************************************************************************
35  !
36  ! variables declaration
37  !     
38  CHARACTER*20,DIMENSION(:),POINTER :: Ncdf_varname => NULL()
39  CHARACTER*20 :: vert_coord_name
40  CHARACTER*1 :: posvar
41  CHARACTER*100 :: Child_file,Childcoordinates,varname,Child_Bathy_Level,Child_Bathy_Meter   
42  REAL*8, POINTER, DIMENSION(:,:) :: lonChild,latChild => NULL()
43  REAL*8, POINTER, DIMENSION(:,:) :: lonParent,latParent => NULL()
44  REAL*8, POINTER, DIMENSION(:,:,:) :: tabvar3d,tabinterp3d,mask => NULL()
45  REAL*8, POINTER, DIMENSION(:,:,:) :: fmask,fse3u,fse3v,fse3t => NULL()
46  REAL*8, POINTER, DIMENSION(:,:,:,:) :: un,ub,vn,vb,tn,tb,sn,sb,e3t_n,e3t_b => NULL()
47  REAL*8, POINTER, DIMENSION(:,:,:) :: sshn,sshb => NULL()     
48  REAL*8, POINTER, DIMENSION(:,:,:,:) :: tabinterp4d,tabvar1,tabvar2,tabvar3 => NULL()
49  REAL*8, POINTER, DIMENSION(:) :: timedepth_temp => NULL()
50  REAL*8, POINTER, DIMENSION(:) :: tabtemp1D,nav_lev => NULL()
51  INTEGER, POINTER, DIMENSION(:) :: tabtemp1DInt => NULL()
52  REAL*8, POINTER, DIMENSION(:,:) :: tabtemp2D,zwf => NULL()
53  REAL*8, POINTER, DIMENSION(:,:) :: e1f,e2f,e1v,e2v,e1u,e2u,e1t,e2t => NULL()
54  REAL*8, POINTER, DIMENSION(:,:,:,:) :: tabtemp4D => NULL()
55  INTEGER,DIMENSION(:),POINTER :: src_add,dst_add  => NULL()
56  REAL*8,DIMENSION(:,:),POINTER :: matrix,bathy_G0 => NULL()
57  LOGICAL,DIMENSION(:,:,:),POINTER :: Tmask => NULL()
58  LOGICAL,DIMENSION(:,:),POINTER :: masksrc => NULL()
59  LOGICAL, DIMENSION(:,:,:), POINTER :: detected_pts
60  LOGICAL :: Interpolation,Extrapolation,Pacifique,op
61  REAL*8 :: za1,za0,zsur,zacr,zkth,zdepth,zdepwp,zmin,zmax,zdiff,ze3tp,ze3wp       
62  INTEGER :: narg,iargc,ncid,x,y,z,z_a,x_a,y_a,z_b,nbvert_lev,m
63  REAL*8 :: now_wght,before_wght
64  INTEGER :: i,j,k,status,ji,jj
65  CHARACTER(len=20),DIMENSION(4) :: dimnames
66  CHARACTER(len=80) :: namelistname
67  TYPE(Coordinates) :: G0,G1
68  INTEGER :: jpi,jpj,jpk,jpni,jpnj,jpnij,jpiglo,jpjglo,nlcit,nlcjt,nldit
69  INTEGER :: nldjt,nleit,nlejt,nimppt,njmppt
70  REAL*8 :: tabtemp0dreal
71  INTEGER :: tabtemp0dint
72  CHARACTER(len=20) :: timedimname
73
74  LOGICAL, PARAMETER :: conservation = .FALSE.
75  !     
76  ! Variables for dimg
77  !
78!  INTEGER :: ino0, it0, ipcg0, isor0, itke0,nfice,nfbulk
79!  INTEGER :: irecl8, irec,ndastp,narea,nsolv
80!  INTEGER :: jk,kt            ! dummy loop indices
81!  INTEGER :: inum             ! temporary logical unit
82!  INTEGER :: ios1 , ios2      ! flag for ice and bulk in the current run
83!  INTEGER :: ios3             ! flag for free surface.  0 = none 1 = yes.  0 = none 1 = yes
84!  INTEGER :: ios4             ! flag for coupled (1) or not (0)   
85  !       
86  narg = iargc()
87  IF (narg == 0) THEN
88     namelistname = 'namelist.input'
89  ELSE
90     CALL getarg(1,namelistname)
91  ENDIF
92  !
93  ! read input file
94  !
95  CALL read_namelist(namelistname)
96  !     
97  IF(TRIM(restart_file) == '') THEN
98     WRITE(*,*) 'no restart file specified in ',TRIM(namelistname)
99     STOP
100  END IF
101
102  IF (iom_activated) THEN
103     timedimname = 'time_counter'
104  ELSE
105     timedimname='time'
106  ENDIF
107
108  !
109  WRITE(*,*) ''
110  WRITE(*,*) 'Interpolation of restart file : ',TRIM(restart_file)
111  WRITE(*,*) ''
112  !
113  CALL Read_Ncdf_VarName(restart_file,Ncdf_varname)
114  !       
115  CALL set_child_name(parent_coordinate_file,Childcoordinates)   
116  IF( TRIM(parent_bathy_level) /= '' )   CALL set_child_name(parent_bathy_level,Child_Bathy_Level) 
117  IF( TRIM(parent_bathy_meter) /= '' )   CALL set_child_name(parent_bathy_meter,Child_Bathy_Meter)   
118  !
119  ! create this file
120  !
121  CALL set_child_name(restart_file,Child_file)
122  status = nf90_create(Child_file,NF90_WRITE,ncid)
123  status = nf90_close(ncid)
124  WRITE(*,*) 'Child grid restart file name = ',TRIM(Child_file)     
125  WRITE(*,*) ''
126
127  !
128  ! read dimensions in parent restart file
129  !
130  CALL Read_Ncdf_dim('x',restart_file,x)
131  CALL Read_Ncdf_dim('y',restart_file,y) 
132  CALL Read_Ncdf_dim('nav_lev',restart_file,z)
133  IF (.NOT.iom_activated) THEN
134     CALL Read_Ncdf_dim('x_a',restart_file,x_a)
135     CALL Read_Ncdf_dim('y_a',restart_file,y_a)
136     CALL Read_Ncdf_dim('z_a',restart_file,z_a)
137     CALL Read_Ncdf_dim('z_b',restart_file,z_b)
138  ENDIF
139
140  IF( z .NE. N ) THEN
141     WRITE(*,*) '***'
142     WRITE(*,*) 'Number of vertical levels doesn t match between namelist and restart file'
143     WRITE(*,*) 'Please check the values in namelist file'
144     STOP
145  ENDIF
146  !
147  ! mask initialization for extrapolation and interpolation
148  !
149  WRITE(*,*) 'mask initialisation on coarse and fine grids'
150  !           
151  status = Read_Local_Coordinates(parent_coordinate_file,G0,(/jpizoom,jpjzoom/),(/x,y/))
152  status = Read_Coordinates(Childcoordinates,G1,Pacifique)
153  !
154  !longitude modification if child domain covers Pacific ocean area
155  !     
156  IF( Pacifique ) THEN
157     !
158     WHERE( G0%nav_lon < 0 )
159        G0%nav_lon = G0%nav_lon + 360.
160     END WHERE
161     !             
162     WHERE( G1%nav_lon < 0 )
163        G1%nav_lon = G1%nav_lon + 360.
164     END WHERE
165     !
166  ENDIF
167  !
168  ! one needs bathy_level
169  IF( TRIM(parent_bathy_level) /= '' ) THEN
170     status = Read_bathy_level(TRIM(parent_bathy_level),G0)
171     status = Read_bathy_level(TRIM(child_bathy_level),G1)
172  ELSE
173     status = read_bathy_meter(TRIM(parent_bathy_meter),G0)
174     status = read_bathy_meter(TRIM(child_bathy_meter),G1)
175     CALL meter_to_levels(G0)
176     CALL meter_to_levels(G1)
177  ENDIF
178  ! get masks
179  CALL Init_mask(parent_bathy_level,G0,x,y)
180  CALL Init_mask(child_bathy_level,G1,1,1)
181
182  G0%tmask = 1.   
183
184  DO k=1,z
185     ALLOCATE(tabvar1(x,y,1,1))
186     CALL Read_Ncdf_var('sn',TRIM(restart_file),tabvar1,1,k)
187     WHERE( tabvar1(:,:,1,1) == 0. ) 
188        G0%tmask(:,:,k) = 0.
189     END WHERE
190     DEALLOCATE(tabvar1)
191  END DO
192  !
193  G0%umask(1:x-1,:,:) = G0%tmask(1:x-1,:,:)*G0%tmask(2:x,:,:)
194  G0%umask(x,:,:)     = G0%tmask(x,:,:)
195  G0%vmask(:,1:y-1,:) = G0%tmask(:,1:y-1,:)*G0%tmask(:,2:y,:)
196  G0%vmask(:,y,:)     = G0%tmask(:,y,:)
197  !     
198  G0%fmask(1:x-1,1:y-1,:) = G0%tmask(1:x-1,1:y-1,:)*G0%tmask(2:x,1:y-1,:)* &
199     &                      G0%tmask(1:x-1,2:y,:)*G0%tmask(2:x,2:y,:) 
200  G0%fmask(x,:,:) = G0%tmask(x,:,:)
201  G0%fmask(:,y,:) = G0%tmask(:,y,:)
202  !
203  !
204  ! write dimensions in output file
205  WRITE(*,*) 'write dimensions'
206  !         
207  CALL Write_Ncdf_dim('x',Child_file,nxfin)
208  CALL Write_Ncdf_dim('y',Child_file,nyfin)
209  CALL Write_Ncdf_dim('nav_lev',Child_file,z)
210  CALL Write_Ncdf_dim(TRIM(timedimname),Child_file,0) 
211  IF (.NOT.iom_activated) THEN
212     CALL Write_Ncdf_dim('x_a',Child_file,x_a)
213     CALL Write_Ncdf_dim('y_a',Child_file,y_a)
214     CALL Write_Ncdf_dim('z_a',Child_file,z_a)
215     CALL Write_Ncdf_dim('z_b',Child_file,z_b)
216  ENDIF
217  !
218  !
219  !
220  !
221  DO ji = 1,SIZE(Ncdf_varname)     
222     !     
223     ! loop on variables names
224     varname = TRIM(Ncdf_varname(ji))
225     WRITE(*,*) 'var = ',TRIM(varname)     
226     !     
227     SELECT CASE (TRIM(varname))
228        !
229     CASE('nav_lon')
230        CALL Read_Ncdf_var('nav_lon',TRIM(Childcoordinates),tabtemp2D) 
231        CALL Write_Ncdf_var('nav_lon',(/'x','y'/),Child_file,tabtemp2D,'float')
232        CALL Copy_Ncdf_att('nav_lon',TRIM(restart_file),Child_file,MINVAL(tabtemp2D),MAXVAL(tabtemp2D))
233        DEALLOCATE(tabtemp2D)
234        Interpolation = .FALSE.
235        !       
236     CASE('nav_lat')             
237        CALL Read_Ncdf_var('nav_lat',TRIM(Childcoordinates),tabtemp2D) 
238        CALL Write_Ncdf_var('nav_lat',(/'x','y'/),Child_file,tabtemp2D,'float')
239        CALL Copy_Ncdf_att('nav_lat',TRIM(restart_file),Child_file,MINVAL(tabtemp2D),MAXVAL(tabtemp2D)) 
240        DEALLOCATE(tabtemp2D)
241        Interpolation = .FALSE.
242        !
243     CASE('nav_lev')
244        CALL Read_Ncdf_var('nav_lev',TRIM(restart_file),nav_lev) 
245        CALL Write_Ncdf_var('nav_lev','nav_lev',Child_file,nav_lev,'float')
246        CALL Copy_Ncdf_att('nav_lev',TRIM(restart_file),Child_file)     
247        Interpolation = .FALSE.
248        !
249     CASE('time')
250        CALL Read_Ncdf_var('time',TRIM(restart_file),tabtemp1D) 
251        CALL Write_Ncdf_var('time',TRIM(timedimname),Child_file,tabtemp1D,'float')
252        CALL Copy_Ncdf_att('time',TRIM(restart_file),Child_file) 
253        DEALLOCATE(tabtemp1D)
254        Interpolation = .FALSE.
255        !
256     CASE('time_counter')
257        CALL Read_Ncdf_var('time_counter',TRIM(restart_file),tabtemp1D) 
258        tabtemp1D = tabtemp1D * rhot
259        CALL Write_Ncdf_var('time_counter',TRIM(timedimname),Child_file,tabtemp1D,'double')
260        CALL Copy_Ncdf_att('time_counter',TRIM(restart_file),Child_file) 
261        DEALLOCATE(tabtemp1D)
262        Interpolation = .FALSE.
263        !
264     CASE('kt','ndastp','adatrj','ntime','nn_fsbc','rdt') 
265        IF (iom_activated) THEN
266           CALL Read_Ncdf_var(TRIM(varname),TRIM(restart_file),tabtemp0dreal) 
267           SELECT CASE (TRIM(varname))
268           CASE('rdt')
269              tabtemp0dreal = tabtemp0dreal / rhot
270           CASE('kt')
271              tabtemp0dreal = tabtemp0dreal * rhot
272           END SELECT
273           CALL Write_Ncdf_var(TRIM(varname),Child_file,tabtemp0dreal,'double')
274        ELSE
275           CALL Read_Ncdf_var(TRIM(varname),TRIM(restart_file),tabtemp4D) 
276           dimnames(1)='x_a'
277           dimnames(2)='y_a'
278           dimnames(3)='z_b'
279           dimnames(4)=TRIM(timedimname)           
280           CALL Write_Ncdf_var(TRIM(varname),dimnames,Child_file,tabtemp4D,'double')
281           DEALLOCATE(tabtemp4D)
282        ENDIF
283        CALL Copy_Ncdf_att(TRIM(varname),TRIM(restart_file),Child_file) 
284        Interpolation = .FALSE.
285        !
286     CASE('frc_v','frc_t','frc_s') 
287        CALL Read_Ncdf_var(varname,TRIM(restart_ice_file),tabtemp0dreal) 
288        CALL Write_Ncdf_var(varname,Child_file,tabtemp0dreal,'double')
289        CALL Copy_Ncdf_att(varname,TRIM(restart_ice_file),Child_file) 
290        Interpolation = .FALSE.
291        !
292        ! Variable interpolation according to their position on grid
293        !                                   
294     CASE('ssu_m','utau_b','un_bf','un','ub') 
295        IF( Get_NbDims(TRIM(varname),TRIM(restart_file)) == 4 ) THEN
296           vert_coord_name = 'nav_lev'
297        ELSEIF( Get_NbDims(TRIM(varname),TRIM(restart_file)) == 3 ) THEN
298           vert_coord_name = '1'
299        ENDIF
300        posvar='U'
301        Interpolation = .TRUE. 
302        !
303     CASE('ssv_m','vtau_b','vn_bf','vn','vb') 
304        IF( Get_NbDims(TRIM(varname),TRIM(restart_file)) == 4 ) THEN
305           vert_coord_name = 'nav_lev'
306        ELSEIF( Get_NbDims(TRIM(varname),TRIM(restart_file)) == 3 ) THEN
307           vert_coord_name = '1'
308        ENDIF
309        posvar='V'
310        Interpolation = .TRUE. 
311        !
312     CASE DEFAULT
313        IF( Get_NbDims(TRIM(varname),TRIM(restart_file)) == 4 ) THEN
314           vert_coord_name = 'nav_lev'
315        ELSEIF( Get_NbDims(TRIM(varname),TRIM(restart_file)) == 3 ) THEN
316           vert_coord_name = '1'
317        ENDIF
318        posvar='T'
319        Interpolation = .TRUE.           
320        !     
321     END SELECT
322     !     
323     ! --- start interpolation --- !
324     IF( Interpolation ) THEN
325        !       
326        IF( vert_coord_name == '1' ) THEN
327           nbvert_lev = 1
328        ELSE
329           nbvert_lev = z
330        ENDIF
331        !
332
333        ALLOCATE(detected_pts(SIZE(G0%tmask,1),SIZE(G0%tmask,2),nbvert_lev))                             
334        ALLOCATE(tabvar1(x,y,1,2))
335        ALLOCATE(tabvar2(x,y,1,1))
336        ALLOCATE(tabvar3(x,y,1,1))
337        ALLOCATE(masksrc(x,y))
338        ALLOCATE(tabinterp4d(nxfin,nyfin,1,1)) 
339
340        !       
341        DO n = 1,nbvert_lev
342           !
343           WRITE(*,*) 'interpolate/extrapolate for vertical level = ',n   
344           !                           
345           CALL Read_Ncdf_var(varname,TRIM(restart_file),tabvar1,1,n)
346           IF(n==1) THEN
347              !                           
348           ELSE IF (n==2) THEN
349              tabvar2(:,:,:,1) = tabvar1(:,:,:,2)
350           ELSE
351              tabvar3(:,:,:,1) = tabvar2(:,:,:,1)
352              tabvar2(:,:,:,1) = tabvar1(:,:,:,2)
353
354           ENDIF
355           !                           
356           SELECT CASE(posvar)
357              !
358           CASE('T')
359              !
360              IF(MAXVAL(G1%tmask(:,:,n)) == 0.) THEN           
361                 tabinterp4d = 0.0
362                 WRITE(*,*) 'only land points on level ',n                     
363              ELSE                       
364
365                 CALL extrap_detect(G0,G1,detected_pts(:,:,n),n)                                           
366
367                 CALL correct_field(detected_pts(:,:,n),tabvar1,tabvar2,tabvar3,G0,nav_lev,masksrc,n)                                 
368
369                 ! for the following variables, you do not want to mask the values
370                 IF(  TRIM(varname) == 'e3t_n' .OR. TRIM(varname) == 'e3t_b' .OR. TRIM(varname) == 'e3t_m' .OR. &
371                    & TRIM(varname) == 'fraqsr_1lev' .OR. TRIM(varname) == 'frq_m' .OR. TRIM(varname) == 'surf_ini' ) THEN
372                    masksrc(:,:) = .TRUE.
373                 ENDIF
374                 
375                 SELECT CASE(TRIM(interp_type))
376                 CASE('bilinear')                                                       
377                    CALL get_remap_matrix(G0%nav_lat,G1%nav_lat, &
378                       G0%nav_lon,G1%nav_lon,masksrc,matrix,src_add,dst_add)
379                    CALL make_remap(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1),nxfin,nyfin, &
380                       matrix,src_add,dst_add)     
381                 CASE('bicubic')                                   
382                    CALL get_remap_bicub(G0%nav_lat,G1%nav_lat, &
383                       G0%nav_lon,G1%nav_lon,masksrc,matrix,src_add,dst_add)
384                    CALL make_bicubic_remap(tabvar1(:,:,1,1),masksrc,tabinterp4d(:,:,1,1),&
385                       nxfin,nyfin,matrix,src_add,dst_add) 
386                 END SELECT
387                 !
388                 IF( conservation ) THEN ! clem: it currently does not work
389                    CALL Correctforconservation(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1), &
390                       G0%e1t,G0%e2t,G1%e1t,G1%e2t,nxfin,nyfin,posvar,imin-jpizoom+1,jmin-jpjzoom+1)
391                 ENDIF
392
393              ENDIF
394             
395              IF( ALL(masksrc) ) THEN
396                 tabinterp4d(:,:,1,1) =  tabinterp4d(:,:,1,1)
397              ELSE
398                 tabinterp4d(:,:,1,1) =  tabinterp4d(:,:,1,1) * G1%tmask(:,:,n)
399              ENDIF
400              !
401           CASE('U')
402              !
403              IF(MAXVAL(G1%umask(:,:,n)) == 0) THEN
404                 tabinterp4d = 0.0
405                 WRITE(*,*) 'only land points on level ',n 
406              ELSE     
407                 !
408                 CALL extrap_detect(G0,G1,detected_pts(:,:,n),n,'U') 
409                 CALL correct_field(detected_pts(:,:,n),tabvar1,tabvar2,tabvar3,G0,nav_lev,masksrc,n,'U')
410                 !                                                           
411                 SELECT CASE(TRIM(interp_type))
412                 CASE('bilinear')                                                       
413                    CALL get_remap_matrix(G0%gphiu,G1%gphiu,   &
414                       G0%glamu,G1%glamu,masksrc,matrix,src_add,dst_add)
415                    CALL make_remap(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1),nxfin,nyfin, &
416                       matrix,src_add,dst_add)     
417                 CASE('bicubic')                                   
418                    CALL get_remap_bicub(G0%gphiu,G1%gphiu,   &
419                       G0%glamu,G1%glamu,masksrc,matrix,src_add,dst_add)
420                    CALL make_bicubic_remap(tabvar1(:,:,1,1),masksrc,tabinterp4d(:,:,1,1),&
421                       nxfin,nyfin,matrix,src_add,dst_add)                       
422                 END SELECT
423                 !                     
424                 IF( conservation ) THEN ! clem: not coded for U
425                    CALL Correctforconservation(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1), &
426                       G0%e1u,G0%e2u,G1%e1u,G1%e2u,nxfin,nyfin,posvar,imin-jpizoom+1,jmin-jpjzoom+1)
427                 ENDIF
428              ENDIF
429
430              tabinterp4d(:,:,1,1) =  tabinterp4d(:,:,1,1) * G1%umask(:,:,n)
431              !
432           CASE('V')
433              !     
434              IF(MAXVAL(G1%vmask(:,:,n)) == 0) THEN
435                 tabinterp4d = 0.0
436                 WRITE(*,*) 'only land points on level ',n 
437              ELSE     
438                 !
439
440                 CALL extrap_detect(G0,G1,detected_pts(:,:,n),n,'V')
441
442                 CALL correct_field(detected_pts(:,:,n),tabvar1,tabvar2,tabvar3,G0,nav_lev,masksrc,n,'V')
443                 !                                                           
444                 SELECT CASE(TRIM(interp_type))
445                 CASE('bilinear')                                                       
446                    CALL get_remap_matrix(G0%gphiv,G1%gphiv,   &
447                       G0%glamv,G1%glamv,masksrc,matrix,src_add,dst_add)
448                    CALL make_remap(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1),nxfin,nyfin, &
449                       matrix,src_add,dst_add)     
450                 CASE('bicubic')                                   
451                    CALL get_remap_bicub(G0%gphiv,G1%gphiv,   &
452                       G0%glamv,G1%glamv,masksrc,matrix,src_add,dst_add)
453                    CALL make_bicubic_remap(tabvar1(:,:,1,1),masksrc,tabinterp4d(:,:,1,1),&
454                       nxfin,nyfin,matrix,src_add,dst_add)                       
455                 END SELECT
456                 !                     
457                 IF( conservation ) THEN ! clem: not coded for V
458                    CALL Correctforconservation(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1), &
459                       G0%e1v,G0%e2v,G1%e1v,G1%e2v,nxfin,nyfin,posvar,imin-jpizoom+1,jmin-jpjzoom+1)
460                 ENDIF
461              ENDIF
462
463              tabinterp4d(:,:,1,1) =  tabinterp4d(:,:,1,1) * G1%vmask(:,:,n)
464              !
465           END SELECT
466           !
467           !
468           dimnames(1)='x'
469           dimnames(2)='y'
470           IF( vert_coord_name == '1' ) THEN
471              dimnames(3)=TRIM(timedimname)
472             
473              ALLOCATE(tabvar3d(SIZE(tabinterp4d,1),SIZE(tabinterp4d,2),SIZE(tabinterp4d,3)))
474              tabvar3d=tabinterp4d(:,:,:,1)
475              CALL Write_Ncdf_var(TRIM(varname),dimnames,Child_file,tabvar3d,1,'double')
476              DEALLOCATE(tabvar3d)
477           ELSE
478              dimnames(3)=vert_coord_name
479              dimnames(4)=TRIM(timedimname)
480             
481              CALL Write_Ncdf_var(TRIM(varname),dimnames,Child_file,tabinterp4d,1,n,'double')
482           ENDIF
483           !             
484           !
485           CALL Copy_Ncdf_att(TRIM(varname),TRIM(restart_file),Child_file)
486           !
487           !
488           IF(ASSOCIATED(matrix)) DEALLOCATE(matrix,src_add,dst_add)   
489           !
490        END DO
491        !
492        DEALLOCATE(detected_pts)
493        DEALLOCATE(tabinterp4d)
494        DEALLOCATE(tabvar1,tabvar2,tabvar3)
495        DEALLOCATE(masksrc) 
496        !                     
497     ENDIF
498
499  END DO
500
501  ! change the before fields
502  IF(rhot == 1) THEN
503     WRITE(*,*) ''   
504     WRITE(*,*) 'no time interpolation (time refinement ratio = 1)'
505  ELSE   
506     now_wght = (rhot-1.)/rhot
507     before_wght = 1./rhot
508     !   
509     ! --- 4D variables --- !               
510     CALL Read_Ncdf_var('un',Child_file,un)     
511     CALL Read_Ncdf_var('vn',Child_file,vn)
512     CALL Read_Ncdf_var('ub',Child_file,ub)
513     CALL Read_Ncdf_var('vb',Child_file,vb)
514     ub = now_wght*un + before_wght*ub
515     vb = now_wght*vn + before_wght*vb
516     !
517     CALL Read_Ncdf_var('tb',Child_file,tb)
518     CALL Read_Ncdf_var('tn',Child_file,tn)
519     tb = now_wght*tn + before_wght*tb
520     !
521     CALL Read_Ncdf_var('sb',Child_file,sb)
522     CALL Read_Ncdf_var('sn',Child_file,sn)           
523     sb = now_wght*sn + before_wght*sb
524     !
525     CALL Read_Ncdf_var('e3t_b',Child_file,e3t_b)
526     CALL Read_Ncdf_var('e3t_n',Child_file,e3t_n)           
527     e3t_b = now_wght*e3t_n + before_wght*e3t_b
528     !
529     dimnames(1)='x'
530     dimnames(2)='y'
531     dimnames(3)='nav_lev'
532     dimnames(4)=TRIM(timedimname)
533     CALL Write_Ncdf_var('ub',dimnames,Child_file,ub,'double')
534     CALL Write_Ncdf_var('vb',dimnames,Child_file,vb,'double')
535     CALL Write_Ncdf_var('tb',dimnames,Child_file,tb,'double')
536     CALL Write_Ncdf_var('sb',dimnames,Child_file,sb,'double')
537     CALL Write_Ncdf_var('e3t_b',dimnames,Child_file,e3t_b,'double')
538     !
539     DEALLOCATE(un,ub,vn,vb,tn,tb,sn,sb,e3t_n,e3t_b)
540     !----------------------               
541     !
542     ! --- 3D variables --- !       
543     CALL Read_Ncdf_var('sshb',Child_file,sshb)
544     CALL Read_Ncdf_var('sshn',Child_file,sshn)
545     sshb = now_wght*sshn + before_wght*sshb
546     !
547     dimnames(1)='x'
548     dimnames(2)='y'
549     dimnames(3)=TRIM(timedimname)
550     CALL Write_Ncdf_var('sshb',dimnames,Child_file,sshb,'double')
551     !               
552     DEALLOCATE(sshn,sshb) 
553     !----------------------               
554     !   
555  ENDIF
556  !
557  !
558  WRITE(*,*) ' '
559  WRITE(*,*) ' --- list of all the variables that have been interpolated --- '
560  WRITE(*,*) Ncdf_varname
561  WRITE(*,*) ' '
562  WRITE(*,*) '******* restart file successfully created *******' 
563  WRITE(*,*) ' '
564  !
565  STOP
566END PROGRAM
567
568
Note: See TracBrowser for help on using the repository browser.