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 @ 12253

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

BGC restart should now work in the nesting tool. All of that must of course be tested in realistic simulations. Waiting for the feedbacks...

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