[2136] | 1 | ! |
---|
| 2 | !************************************************************************ |
---|
| 3 | ! Fortran 95 OPA Nesting tools * |
---|
| 4 | ! * |
---|
| 5 | ! Copyright (C) 2005 Florian Lemari�(Florian.Lemarie@imag.fr) * |
---|
| 6 | ! * |
---|
| 7 | !************************************************************************ |
---|
| 8 | ! |
---|
| 9 | PROGRAM 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 | ! |
---|
| 20 | IMPLICIT NONE |
---|
| 21 | ! |
---|
| 22 | !************************************************************************ |
---|
| 23 | ! * |
---|
| 24 | ! PROGRAM CREATE_RSTRT * |
---|
| 25 | ! * |
---|
| 26 | ! program to interpolate parent grid restart file to child grid * |
---|
| 27 | ! * |
---|
| 28 | ! * |
---|
| 29 | !Interpolation is carried out using bilinear interpolation * |
---|
| 30 | !routine from SCRIP package * |
---|
| 31 | ! * |
---|
| 32 | !http://climate.lanl.gov/Software/SCRIP/ * |
---|
| 33 | !************************************************************************ |
---|
| 34 | ! |
---|
| 35 | ! variables declaration |
---|
| 36 | ! |
---|
| 37 | CHARACTER*20,DIMENSION(:),POINTER :: Ncdf_varname => NULL() |
---|
| 38 | CHARACTER*20 :: vert_coord_name |
---|
| 39 | CHARACTER*1 :: posvar |
---|
| 40 | CHARACTER*100 :: Child_file,Childcoordinates,varname,Childbathy,Childbathymeter |
---|
| 41 | REAL*8, POINTER, DIMENSION(:,:) :: lonChild,latChild => NULL() |
---|
| 42 | REAL*8, POINTER, DIMENSION(:,:) :: lonParent,latParent => NULL() |
---|
| 43 | REAL*8, POINTER, DIMENSION(:,:,:) :: tabvar3d,tabinterp3d,mask => NULL() |
---|
| 44 | REAL*8, POINTER, DIMENSION(:,:,:) :: fmask,fse3u,fse3v,fse3t => NULL() |
---|
| 45 | REAL*8, POINTER, DIMENSION(:,:,:,:) :: un,ub,vn,vb,tn,tb => NULL() |
---|
| 46 | REAL*8, POINTER, DIMENSION(:,:,:,:) :: sshn,sshb,sb,sn,gcx,gcxb => NULL() |
---|
| 47 | REAL*8, POINTER, DIMENSION(:,:,:,:) :: tabinterp4d,tabvar1,tabvar2,tabvar3 => NULL() |
---|
| 48 | REAL*8, POINTER, DIMENSION(:,:,:,:) :: rotn,rotb,hdivn,hdivb => 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,t,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 | ! |
---|
| 75 | ! Variables for dimg |
---|
| 76 | ! |
---|
| 77 | INTEGER :: ino0, it0, ipcg0, isor0, itke0,nfice,nfbulk |
---|
| 78 | INTEGER :: irecl8, irec,ndastp,narea,nsolv |
---|
| 79 | INTEGER :: jk,kt ! dummy loop indices |
---|
| 80 | INTEGER :: inum ! temporary logical unit |
---|
| 81 | INTEGER :: ios1 , ios2 ! flag for ice and bulk in the current run |
---|
| 82 | INTEGER :: ios3 ! flag for free surface. 0 = none 1 = yes. 0 = none 1 = yes |
---|
| 83 | INTEGER :: ios4 ! flag for coupled (1) or not (0) |
---|
| 84 | ! |
---|
| 85 | narg = iargc() |
---|
| 86 | IF (narg == 0) THEN |
---|
| 87 | namelistname = 'namelist.input' |
---|
| 88 | ELSE |
---|
| 89 | CALL getarg(1,namelistname) |
---|
| 90 | ENDIF |
---|
| 91 | ! |
---|
| 92 | ! read input file |
---|
| 93 | ! |
---|
| 94 | CALL read_namelist(namelistname) |
---|
| 95 | ! |
---|
| 96 | IF(TRIM(restart_file) == '/NULL') THEN |
---|
| 97 | WRITE(*,*) 'no restart file specified in ',TRIM(namelistname) |
---|
| 98 | STOP |
---|
| 99 | END IF |
---|
| 100 | |
---|
| 101 | IF (iom_activated) THEN |
---|
| 102 | timedimname = 't' |
---|
| 103 | ELSE |
---|
| 104 | timedimname='time' |
---|
| 105 | ENDIF |
---|
| 106 | |
---|
| 107 | ! |
---|
| 108 | WRITE(*,*) '' |
---|
| 109 | WRITE(*,*) 'Interpolation of restart file : ',TRIM(restart_file) |
---|
| 110 | WRITE(*,*) '' |
---|
| 111 | ! |
---|
| 112 | CALL Read_Ncdf_VarName(restart_file,Ncdf_varname) |
---|
| 113 | ! |
---|
| 114 | CALL set_child_name(parent_coordinate_file,Childcoordinates) |
---|
[2455] | 115 | CALL set_child_name(parent_meshmask_file,Childbathy) |
---|
[2136] | 116 | CALL set_child_name(parent_bathy_meter,Childbathymeter) |
---|
| 117 | ! |
---|
| 118 | ! create this file |
---|
| 119 | ! |
---|
| 120 | IF( .NOT. dimg ) THEN |
---|
| 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 | ENDIF |
---|
| 127 | |
---|
| 128 | ! |
---|
| 129 | ! read dimensions in parent restart file |
---|
| 130 | ! |
---|
| 131 | CALL Read_Ncdf_dim('x',restart_file,x) |
---|
| 132 | CALL Read_Ncdf_dim('y',restart_file,y) |
---|
| 133 | CALL Read_Ncdf_dim('z',restart_file,z) |
---|
| 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 | |
---|
| 139 | IF( z .NE. N ) THEN |
---|
| 140 | WRITE(*,*) '***' |
---|
| 141 | WRITE(*,*) 'Number of vertical levels doesn t match between namelist and restart file' |
---|
| 142 | WRITE(*,*) 'Please check the values in namelist file' |
---|
| 143 | STOP |
---|
| 144 | ENDIF |
---|
| 145 | ! |
---|
| 146 | ! mask initialization for extrapolation and interpolation |
---|
| 147 | ! |
---|
| 148 | WRITE(*,*) 'mask initialisation on coarse and fine grids' |
---|
| 149 | ! |
---|
| 150 | status = Read_Local_Coordinates(parent_coordinate_file,G0,(/jpizoom,jpjzoom/),(/x,y/)) |
---|
| 151 | status = Read_Coordinates(Childcoordinates,G1,Pacifique) |
---|
| 152 | ! |
---|
| 153 | !longitude modification if child domain covers Pacific ocean area |
---|
| 154 | ! |
---|
| 155 | IF( Pacifique ) THEN |
---|
| 156 | ! |
---|
| 157 | WHERE( G0%nav_lon < 0 ) |
---|
| 158 | G0%nav_lon = G0%nav_lon + 360. |
---|
| 159 | END WHERE |
---|
| 160 | ! |
---|
| 161 | WHERE( G1%nav_lon < 0 ) |
---|
| 162 | G1%nav_lon = G1%nav_lon + 360. |
---|
| 163 | END WHERE |
---|
| 164 | ! |
---|
| 165 | ENDIF |
---|
| 166 | ! |
---|
[2455] | 167 | CALL Init_mask(parent_meshmask_file,G0,x,y) |
---|
[2136] | 168 | CALL Init_mask(childbathy,G1,1,1) |
---|
| 169 | |
---|
| 170 | G0%tmask = 1. |
---|
| 171 | |
---|
| 172 | DO k=1,z |
---|
| 173 | ALLOCATE(tabvar1(x,y,1,1)) |
---|
| 174 | CALL Read_Ncdf_var('sn',TRIM(restart_file),tabvar1,1,k) |
---|
| 175 | WHERE( tabvar1(:,:,1,1) == 0. ) |
---|
| 176 | G0%tmask(:,:,k) = 0. |
---|
| 177 | END WHERE |
---|
| 178 | DEALLOCATE(tabvar1) |
---|
| 179 | END DO |
---|
| 180 | ! |
---|
| 181 | G0%umask(1:x-1,:,:) = G0%tmask(1:x-1,:,:)*G0%tmask(2:x,:,:) |
---|
| 182 | G0%vmask(:,1:y-1,:) = G0%tmask(:,1:y-1,:)*G0%tmask(:,2:y,:) |
---|
| 183 | ! |
---|
| 184 | G0%umask(x,:,:) = G0%tmask(x,:,:) |
---|
| 185 | G0%vmask(:,y,:) = G0%tmask(:,y,:) |
---|
| 186 | ! |
---|
| 187 | G0%fmask(1:x-1,1:y-1,:) = G0%tmask(1:x-1,1:y-1,:)*G0%tmask(2:x,1:y-1,:)* & |
---|
| 188 | G0%tmask(1:x-1,2:y,:)*G0%tmask(2:x,2:y,:) |
---|
| 189 | ! |
---|
| 190 | G0%fmask(x,:,:) = G0%tmask(x,:,:) |
---|
| 191 | G0%fmask(:,y,:) = G0%tmask(:,y,:) |
---|
| 192 | ! |
---|
| 193 | ! ***** *** ** ** ****** |
---|
| 194 | ! * * * * * * * * |
---|
| 195 | ! * * * * * * * *** |
---|
| 196 | ! * * * * * * * |
---|
| 197 | ! ***** *** * * ****** |
---|
| 198 | ! |
---|
| 199 | IF(dimg) THEN |
---|
| 200 | ! |
---|
| 201 | WRITE(*,*) 'create dimg restart file' |
---|
| 202 | DO m = 7,1000 |
---|
| 203 | INQUIRE(Unit=inum,Opened=op) |
---|
| 204 | IF( .NOT. op ) THEN |
---|
| 205 | inum = m |
---|
| 206 | EXIT |
---|
| 207 | ENDIF |
---|
| 208 | ENDDO |
---|
| 209 | ! |
---|
| 210 | inum = 11 |
---|
| 211 | irecl8 = nxfin * nyfin * 8 |
---|
| 212 | OPEN(inum,FILE=TRIM(dimg_output_file),FORM='UNFORMATTED', & |
---|
| 213 | ACCESS='DIRECT',RECL=irecl8) |
---|
| 214 | ! |
---|
| 215 | CALL Read_Ncdf_var('info',TRIM(restart_file),tabtemp4D) |
---|
| 216 | ! |
---|
| 217 | ino0 = NINT(tabtemp4D(1,1,1,1)) |
---|
| 218 | it0 = NINT(tabtemp4D(1,1,2,1))*rhot |
---|
| 219 | WRITE(*,*) 'restart file created for kt = ',it0 |
---|
| 220 | ipcg0 = NINT(tabtemp4D(1,1,3,1)) |
---|
| 221 | isor0 = NINT(tabtemp4D(1,1,4,1)) |
---|
| 222 | itke0 = 0 |
---|
| 223 | IF(tabtemp4D(1,1,5,1)==1.) itke0 = 1 |
---|
| 224 | ndastp = NINT(tabtemp4D(1,1,6,1)) |
---|
| 225 | ! number of elapsed days since the begining of the run |
---|
| 226 | ! |
---|
| 227 | DEALLOCATE(tabtemp4D) |
---|
| 228 | ! |
---|
| 229 | IF (isor0 + 1 == 3) THEN |
---|
| 230 | isor0 = 2 |
---|
| 231 | ipcg0 = 2 |
---|
| 232 | ENDIF |
---|
| 233 | ! |
---|
| 234 | CALL Read_Ncdf_var('nfice',TRIM(restart_file),tabtemp4D) |
---|
| 235 | nfice = NINT(tabtemp4D(1,1,1,1)) |
---|
| 236 | DEALLOCATE(tabtemp4D) |
---|
| 237 | ! |
---|
| 238 | CALL Read_Ncdf_var('nfbulk',TRIM(restart_file),tabtemp4D) |
---|
| 239 | nfbulk = NINT(tabtemp4D(1,1,1,1)) |
---|
| 240 | DEALLOCATE(tabtemp4D) |
---|
| 241 | ! |
---|
| 242 | nfice = 5 |
---|
| 243 | nfbulk = 5 |
---|
| 244 | ! |
---|
| 245 | ios1 = 0 |
---|
| 246 | ios2 = 0 |
---|
| 247 | ios3 = 1 ! flag for free surface. 0 = none 1 = yes. 0 |
---|
| 248 | ios4 = 0 |
---|
| 249 | narea = 1 |
---|
| 250 | jpi = nxfin |
---|
| 251 | jpj = nyfin |
---|
| 252 | jpk = z |
---|
| 253 | jpni = 1 |
---|
| 254 | jpnj = 1 |
---|
| 255 | jpnij = 1 |
---|
| 256 | narea = 1 |
---|
| 257 | jpiglo = nxfin |
---|
| 258 | jpjglo = nyfin |
---|
| 259 | nlcit = nxfin |
---|
| 260 | nlcjt = nyfin |
---|
| 261 | nldit = 1 |
---|
| 262 | nldjt = 1 |
---|
| 263 | nleit = nxfin |
---|
| 264 | nlejt = nyfin |
---|
| 265 | nimppt = 1 |
---|
| 266 | njmppt = 1 |
---|
| 267 | ! |
---|
| 268 | PRINT*,'jpi = ',jpi |
---|
| 269 | PRINT*,'jpj = ',jpj |
---|
| 270 | PRINT*,'jpk = ',jpk |
---|
| 271 | PRINT*,'nfice = ',nfice |
---|
| 272 | PRINT*,'nfbulk = ',nfbulk |
---|
| 273 | PRINT*,'ndastp = ',ndastp |
---|
| 274 | ! |
---|
| 275 | WRITE(inum,REC=1) irecl8, ino0, it0, isor0, ipcg0, itke0, & |
---|
| 276 | nfice, nfbulk , ios1, ios2, ios3, ios4, & |
---|
| 277 | ndastp, adatrj, jpi, jpj, jpk, & |
---|
| 278 | jpni, jpnj, jpnij, narea, jpiglo, jpjglo, & |
---|
| 279 | nlcit, nlcjt, nldit, nldjt, nleit, nlejt, nimppt, njmppt |
---|
| 280 | ! |
---|
| 281 | irec = 2 |
---|
| 282 | ! |
---|
| 283 | ! ***** ***** ****** |
---|
| 284 | ! * * * * |
---|
| 285 | ! * * * *** |
---|
| 286 | ! * * * * |
---|
| 287 | ! ***** ***** * |
---|
| 288 | ! |
---|
| 289 | ELSE |
---|
| 290 | |
---|
| 291 | |
---|
| 292 | ! |
---|
| 293 | ! write dimensions in output file |
---|
| 294 | ! |
---|
| 295 | CALL Write_Ncdf_dim('x',Child_file,nxfin) |
---|
| 296 | CALL Write_Ncdf_dim('y',Child_file,nyfin) |
---|
| 297 | CALL Write_Ncdf_dim('z',Child_file,z) |
---|
| 298 | CALL Write_Ncdf_dim(TRIM(timedimname),Child_file,0) |
---|
| 299 | IF (.NOT.iom_activated) THEN |
---|
| 300 | CALL Write_Ncdf_dim('x_a',Child_file,x_a) |
---|
| 301 | CALL Write_Ncdf_dim('y_a',Child_file,y_a) |
---|
| 302 | CALL Write_Ncdf_dim('z_a',Child_file,z_a) |
---|
| 303 | CALL Write_Ncdf_dim('z_b',Child_file,z_b) |
---|
| 304 | ENDIF |
---|
| 305 | ! |
---|
| 306 | ! |
---|
| 307 | ENDIF |
---|
| 308 | |
---|
| 309 | |
---|
| 310 | |
---|
| 311 | |
---|
| 312 | |
---|
| 313 | |
---|
| 314 | ! |
---|
| 315 | ! |
---|
| 316 | ! |
---|
| 317 | DO i = 1,SIZE(Ncdf_varname) |
---|
| 318 | ! |
---|
| 319 | ! loop on variables names |
---|
| 320 | ! |
---|
| 321 | SELECT CASE (TRIM(Ncdf_varname(i))) |
---|
| 322 | ! |
---|
| 323 | !copy nav_lon from child coordinates to output file |
---|
| 324 | ! |
---|
| 325 | CASE('nav_lon') |
---|
| 326 | IF(.NOT. dimg ) THEN |
---|
| 327 | WRITE(*,*) 'copy nav_lon' |
---|
| 328 | CALL Read_Ncdf_var('nav_lon',TRIM(Childcoordinates),tabtemp2D) |
---|
| 329 | CALL Write_Ncdf_var('nav_lon',(/'x','y'/),Child_file,tabtemp2D,'float') |
---|
| 330 | CALL Copy_Ncdf_att('nav_lon',TRIM(restart_file),Child_file, & |
---|
| 331 | MINVAL(tabtemp2D),MAXVAL(tabtemp2D)) |
---|
| 332 | DEALLOCATE(tabtemp2D) |
---|
| 333 | Interpolation = .FALSE. |
---|
| 334 | ENDIF |
---|
| 335 | ! |
---|
| 336 | !copy nav_lat from child coordinates to output file |
---|
| 337 | ! |
---|
| 338 | CASE('nav_lat') |
---|
| 339 | IF(.NOT. dimg ) THEN |
---|
| 340 | WRITE(*,*) 'copy nav_lat' |
---|
| 341 | CALL Read_Ncdf_var('nav_lat',TRIM(Childcoordinates),tabtemp2D) |
---|
| 342 | CALL Write_Ncdf_var('nav_lat',(/'x','y'/),Child_file,tabtemp2D,'float') |
---|
| 343 | CALL Copy_Ncdf_att('nav_lat',TRIM(restart_file),Child_file, & |
---|
| 344 | MINVAL(tabtemp2D),MAXVAL(tabtemp2D)) |
---|
| 345 | DEALLOCATE(tabtemp2D) |
---|
| 346 | Interpolation = .FALSE. |
---|
| 347 | ENDIF |
---|
| 348 | ! |
---|
| 349 | !copy nav_lev from restart_file to output file |
---|
| 350 | ! |
---|
| 351 | CASE('nav_lev') |
---|
| 352 | |
---|
| 353 | WRITE(*,*) 'copy nav_lev' |
---|
| 354 | CALL Read_Ncdf_var('nav_lev',TRIM(restart_file),nav_lev) |
---|
| 355 | IF(.NOT. dimg ) THEN |
---|
| 356 | CALL Write_Ncdf_var('nav_lev','z',Child_file,nav_lev,'float') |
---|
| 357 | CALL Copy_Ncdf_att('nav_lev',TRIM(restart_file),Child_file) |
---|
| 358 | ENDIF |
---|
| 359 | Interpolation = .FALSE. |
---|
| 360 | ! |
---|
| 361 | !copy time from restart_file to output file |
---|
| 362 | ! |
---|
| 363 | CASE('time') |
---|
| 364 | IF(.NOT. dimg ) THEN |
---|
| 365 | WRITE(*,*) 'copy time' |
---|
| 366 | CALL Read_Ncdf_var('time',TRIM(restart_file),tabtemp1D) |
---|
| 367 | CALL Write_Ncdf_var('time',TRIM(timedimname),Child_file,tabtemp1D,'float') |
---|
| 368 | CALL Copy_Ncdf_att('time',TRIM(restart_file),Child_file) |
---|
| 369 | DEALLOCATE(tabtemp1D) |
---|
| 370 | Interpolation = .FALSE. |
---|
| 371 | ENDIF |
---|
| 372 | !copy time from restart_file to output file |
---|
| 373 | ! |
---|
| 374 | CASE('time_counter') |
---|
| 375 | IF(.NOT. dimg ) THEN |
---|
| 376 | WRITE(*,*) 'copy time_counter' |
---|
| 377 | CALL Read_Ncdf_var('time_counter',TRIM(restart_file),tabtemp1D) |
---|
| 378 | tabtemp1D = tabtemp1D * rhot |
---|
| 379 | CALL Write_Ncdf_var('time_counter',TRIM(timedimname),Child_file,tabtemp1D,'double') |
---|
| 380 | CALL Copy_Ncdf_att('time_counter',TRIM(restart_file),Child_file) |
---|
| 381 | DEALLOCATE(tabtemp1D) |
---|
| 382 | Interpolation = .FALSE. |
---|
| 383 | ENDIF |
---|
| 384 | ! |
---|
| 385 | !copy time_steps from restart_file to output file |
---|
| 386 | ! |
---|
| 387 | CASE('time_steps') |
---|
| 388 | IF(.NOT. dimg ) THEN |
---|
| 389 | WRITE(*,*) 'copy time_steps' |
---|
| 390 | CALL Read_Ncdf_var('time_steps',TRIM(restart_file),tabtemp1DInt) |
---|
| 391 | CALL Write_Ncdf_var('time_steps',TRIM(timedimname),Child_file,tabtemp1DInt) |
---|
| 392 | CALL Copy_Ncdf_att('time_steps',TRIM(restart_file),Child_file) |
---|
| 393 | DEALLOCATE(tabtemp1DInt) |
---|
| 394 | Interpolation = .FALSE. |
---|
| 395 | ENDIF |
---|
| 396 | ! |
---|
| 397 | !copy info from restart_file to output file |
---|
| 398 | ! |
---|
| 399 | CASE('info') |
---|
| 400 | IF(.NOT. dimg ) THEN |
---|
| 401 | WRITE(*,*) 'copy info' |
---|
| 402 | CALL Read_Ncdf_var('info',TRIM(restart_file),tabtemp4D) |
---|
| 403 | dimnames(1)='x_a' |
---|
| 404 | dimnames(2)='y_a' |
---|
| 405 | dimnames(3)='z_a' |
---|
| 406 | dimnames(4)=TRIM(timedimname) |
---|
| 407 | CALL Write_Ncdf_var('info',dimnames,Child_file,tabtemp4D,'double') |
---|
| 408 | CALL Copy_Ncdf_att('info',TRIM(restart_file),Child_file) |
---|
| 409 | DEALLOCATE(tabtemp4D) |
---|
| 410 | Interpolation = .FALSE. |
---|
| 411 | ENDIF |
---|
| 412 | ! |
---|
| 413 | CASE('nfice','nfbulk','kt','ndastp','adatrj','rdt','rdttra1') |
---|
| 414 | IF(.NOT. dimg ) THEN |
---|
| 415 | WRITE(*,*) 'copy ',TRIM(Ncdf_varname(i)) |
---|
| 416 | IF (iom_activated) THEN |
---|
| 417 | CALL Read_Ncdf_var(TRIM(Ncdf_varname(i)),TRIM(restart_file),tabtemp0dreal) |
---|
| 418 | SELECT CASE (TRIM(Ncdf_varname(i))) |
---|
| 419 | CASE('rdt','rdttra1') |
---|
| 420 | tabtemp0dreal = tabtemp0dreal/rhot |
---|
| 421 | CASE('kt') |
---|
| 422 | tabtemp0dreal = tabtemp0dreal * rhot |
---|
| 423 | END SELECT |
---|
| 424 | CALL Write_Ncdf_var(TRIM(Ncdf_varname(i)),Child_file,tabtemp0dreal,'double') |
---|
| 425 | ELSE |
---|
| 426 | CALL Read_Ncdf_var(TRIM(Ncdf_varname(i)),TRIM(restart_file),tabtemp4D) |
---|
| 427 | dimnames(1)='x_a' |
---|
| 428 | dimnames(2)='y_a' |
---|
| 429 | dimnames(3)='z_b' |
---|
| 430 | dimnames(4)=TRIM(timedimname) |
---|
| 431 | CALL Write_Ncdf_var(TRIM(Ncdf_varname(i)),dimnames,Child_file,tabtemp4D,'double') |
---|
| 432 | DEALLOCATE(tabtemp4D) |
---|
| 433 | ENDIF |
---|
| 434 | CALL Copy_Ncdf_att(TRIM(Ncdf_varname(i)),TRIM(restart_file),Child_file) |
---|
| 435 | Interpolation = .FALSE. |
---|
| 436 | ENDIF |
---|
| 437 | ! |
---|
| 438 | ! Variable interpolation according to their position on grid |
---|
| 439 | ! |
---|
| 440 | CASE('un','ub') |
---|
| 441 | varname = Ncdf_varname(i) |
---|
| 442 | |
---|
| 443 | IF(TRIM(varname)=='un') irec = 6 * z + 2 |
---|
| 444 | IF(TRIM(varname)=='ub') irec = 2 |
---|
| 445 | |
---|
| 446 | WRITE(*,*) TRIM(varname),'interpolation ...' |
---|
| 447 | vert_coord_name = 'z' |
---|
| 448 | posvar='U' |
---|
| 449 | Interpolation = .TRUE. |
---|
| 450 | ! |
---|
| 451 | CASE('u_io') |
---|
| 452 | varname = Ncdf_varname(i) |
---|
| 453 | WRITE(*,*) TRIM(varname),'interpolation ...' |
---|
| 454 | vert_coord_name = 'z_b' |
---|
| 455 | posvar='U' |
---|
| 456 | Interpolation = .TRUE. |
---|
| 457 | ! |
---|
| 458 | CASE('vn','vb') |
---|
| 459 | varname = Ncdf_varname(i) |
---|
| 460 | |
---|
| 461 | IF(TRIM(varname)=='vn') irec = 7 * z + 2 |
---|
| 462 | IF(TRIM(varname)=='vb') irec = 2 + z |
---|
| 463 | |
---|
| 464 | WRITE(*,*) TRIM(varname),'interpolation ...' |
---|
| 465 | vert_coord_name = 'z' |
---|
| 466 | posvar='V' |
---|
| 467 | Interpolation = .TRUE. |
---|
| 468 | ! |
---|
| 469 | CASE('v_io') |
---|
| 470 | varname = Ncdf_varname(i) |
---|
| 471 | WRITE(*,*) TRIM(varname),'interpolation ...' |
---|
| 472 | vert_coord_name = 'z_b' |
---|
| 473 | posvar='V' |
---|
| 474 | Interpolation = .TRUE. |
---|
| 475 | ! |
---|
| 476 | CASE('gcx','gcxb','sshb','sshn','sst_io','sss_io','gsst') |
---|
| 477 | varname = Ncdf_varname(i) |
---|
| 478 | |
---|
| 479 | IF(TRIM(varname)=='gcx') irec = 12 * z + 2 |
---|
| 480 | IF(TRIM(varname)=='gcxb') irec = 12 * z + 3 |
---|
| 481 | IF(TRIM(varname)=='sshb') irec = 12 * z + 4 |
---|
| 482 | IF(TRIM(varname)=='sshn') irec = 12 * z + 5 |
---|
| 483 | |
---|
| 484 | WRITE(*,*) TRIM(varname),'interpolation ...' |
---|
| 485 | vert_coord_name = 'z_b' |
---|
| 486 | posvar='T' |
---|
| 487 | Interpolation = .TRUE. |
---|
| 488 | |
---|
| 489 | ! |
---|
| 490 | CASE ('tb','sb','sn','tn') |
---|
| 491 | varname = Ncdf_varname(i) |
---|
| 492 | |
---|
| 493 | IF(TRIM(varname)=='sn') irec = 9 * z + 2 |
---|
| 494 | IF(TRIM(varname)=='tn') irec = 8 * z + 2 |
---|
| 495 | IF(TRIM(varname)=='sb') irec = 3 * z + 2 |
---|
| 496 | IF(TRIM(varname)=='tb') irec = 2 * z + 2 |
---|
| 497 | |
---|
| 498 | WRITE(*,*) TRIM(varname),'interpolation ...' |
---|
| 499 | vert_coord_name = 'z' |
---|
| 500 | posvar='T' |
---|
| 501 | Interpolation = .TRUE. |
---|
| 502 | |
---|
| 503 | CASE('en') |
---|
| 504 | varname = Ncdf_varname(i) |
---|
| 505 | irec = 12 * z + 6 |
---|
| 506 | WRITE(*,*) TRIM(varname),'interpolation ...' |
---|
| 507 | vert_coord_name = 'z' |
---|
| 508 | posvar='T' |
---|
| 509 | Interpolation = .TRUE. |
---|
| 510 | |
---|
| 511 | CASE ('rotb','rotn','hdivb','hdivn') |
---|
| 512 | Interpolation = .FALSE. |
---|
| 513 | ! |
---|
| 514 | END SELECT |
---|
| 515 | ! |
---|
| 516 | IF( Interpolation ) THEN |
---|
| 517 | ! |
---|
| 518 | IF( vert_coord_name == 'z') THEN |
---|
| 519 | nbvert_lev = z |
---|
| 520 | ELSE IF( vert_coord_name == 'z_b') THEN |
---|
| 521 | IF (iom_activated) THEN |
---|
| 522 | nbvert_lev=1 |
---|
| 523 | ELSE |
---|
| 524 | nbvert_lev = z_b |
---|
| 525 | ENDIF |
---|
| 526 | END IF |
---|
| 527 | ! |
---|
| 528 | |
---|
| 529 | ! |
---|
| 530 | t = 1 |
---|
| 531 | |
---|
| 532 | ALLOCATE(detected_pts(SIZE(G0%tmask,1),SIZE(G0%tmask,2),nbvert_lev)) |
---|
| 533 | ALLOCATE(tabvar1(x,y,1,2)) |
---|
| 534 | ALLOCATE(tabvar2(x,y,1,1)) |
---|
| 535 | ALLOCATE(tabvar3(x,y,1,1)) |
---|
| 536 | ALLOCATE(masksrc(x,y)) |
---|
| 537 | ALLOCATE(tabinterp4d(nxfin,nyfin,1,1)) |
---|
| 538 | |
---|
| 539 | ! |
---|
| 540 | DO n = 1,nbvert_lev |
---|
| 541 | ! |
---|
| 542 | WRITE(*,*) 'interpolate/extrapolate ', & |
---|
| 543 | TRIM(varname),' for vertical level = ',n |
---|
| 544 | ! |
---|
| 545 | ! |
---|
| 546 | ! If(n==1) then |
---|
| 547 | ! Call Read_Ncdf_var(varname,TRIM(restart_file),tabvar1,t,n) |
---|
| 548 | ! else if (n==2) then |
---|
| 549 | ! Call Read_Ncdf_var(varname,TRIM(restart_file),tabvar2,t,n-1) |
---|
| 550 | ! Call Read_Ncdf_var(varname,TRIM(restart_file),tabvar1,t,n) |
---|
| 551 | ! else |
---|
| 552 | ! Call Read_Ncdf_var(varname,TRIM(restart_file),tabvar3,t,n-2) |
---|
| 553 | ! Call Read_Ncdf_var(varname,TRIM(restart_file),tabvar2,t,n-1) |
---|
| 554 | ! Call Read_Ncdf_var(varname,TRIM(restart_file),tabvar1,t,n) |
---|
| 555 | ! endif |
---|
| 556 | ! |
---|
| 557 | CALL Read_Ncdf_var(varname,TRIM(restart_file),tabvar1,t,n) |
---|
| 558 | IF(n==1) THEN |
---|
| 559 | ! |
---|
| 560 | ELSE IF (n==2) THEN |
---|
| 561 | tabvar2(:,:,:,1) = tabvar1(:,:,:,2) |
---|
| 562 | ELSE |
---|
| 563 | tabvar3(:,:,:,1) = tabvar2(:,:,:,1) |
---|
| 564 | tabvar2(:,:,:,1) = tabvar1(:,:,:,2) |
---|
| 565 | |
---|
| 566 | ENDIF |
---|
| 567 | ! |
---|
| 568 | SELECT CASE(posvar) |
---|
| 569 | ! |
---|
| 570 | CASE('T') |
---|
| 571 | ! |
---|
| 572 | IF(MAXVAL(G1%tmask(:,:,n)) == 0.) THEN |
---|
| 573 | tabinterp4d = 0.0 |
---|
| 574 | WRITE(*,*) 'only land points on level ',n |
---|
| 575 | ELSE |
---|
| 576 | |
---|
| 577 | CALL extrap_detect(G0,G1,detected_pts(:,:,n),n) |
---|
| 578 | |
---|
| 579 | CALL correct_field(detected_pts(:,:,n),tabvar1,tabvar2,& |
---|
| 580 | tabvar3,G0,nav_lev,masksrc,n) |
---|
| 581 | |
---|
| 582 | SELECT CASE(TRIM(interp_type)) |
---|
| 583 | CASE('bilinear') |
---|
| 584 | CALL get_remap_matrix(G0%nav_lat,G1%nav_lat, & |
---|
| 585 | G0%nav_lon,G1%nav_lon,masksrc,matrix,src_add,dst_add) |
---|
| 586 | CALL make_remap(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1),nxfin,nyfin, & |
---|
| 587 | matrix,src_add,dst_add) |
---|
| 588 | CASE('bicubic') |
---|
| 589 | CALL get_remap_bicub(G0%nav_lat,G1%nav_lat, & |
---|
| 590 | G0%nav_lon,G1%nav_lon,masksrc,matrix,src_add,dst_add) |
---|
| 591 | CALL make_bicubic_remap(tabvar1(:,:,1,1),masksrc,tabinterp4d(:,:,1,1),& |
---|
| 592 | nxfin,nyfin,matrix,src_add,dst_add) |
---|
| 593 | END SELECT |
---|
| 594 | ! |
---|
| 595 | CALL Correctforconservation(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1), & |
---|
| 596 | G0%e1t,G0%e2t,G1%e1t,G1%e2t,nxfin,nyfin,posvar,imin-jpizoom+1,jmin-jpjzoom+1) |
---|
| 597 | |
---|
| 598 | |
---|
| 599 | ENDIF |
---|
| 600 | |
---|
| 601 | tabinterp4d(:,:,1,1) = tabinterp4d(:,:,1,1) * G1%tmask(:,:,n) |
---|
| 602 | ! |
---|
| 603 | CASE('U') |
---|
| 604 | ! |
---|
| 605 | IF(MAXVAL(G1%umask(:,:,n)) == 0) THEN |
---|
| 606 | tabinterp4d = 0.0 |
---|
| 607 | WRITE(*,*) 'only land points on level ',n |
---|
| 608 | ELSE |
---|
| 609 | ! |
---|
| 610 | CALL extrap_detect(G0,G1,detected_pts(:,:,n),n,'U') |
---|
| 611 | CALL correct_field(detected_pts(:,:,n),tabvar1,tabvar2,& |
---|
| 612 | tabvar3,G0,nav_lev,masksrc,n,'U') |
---|
| 613 | ! |
---|
| 614 | SELECT CASE(TRIM(interp_type)) |
---|
| 615 | CASE('bilinear') |
---|
| 616 | CALL get_remap_matrix(G0%gphiu,G1%gphiu, & |
---|
| 617 | G0%glamu,G1%glamu,masksrc,matrix,src_add,dst_add) |
---|
| 618 | CALL make_remap(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1),nxfin,nyfin, & |
---|
| 619 | matrix,src_add,dst_add) |
---|
| 620 | CASE('bicubic') |
---|
| 621 | CALL get_remap_bicub(G0%gphiu,G1%gphiu, & |
---|
| 622 | G0%glamu,G1%glamu,masksrc,matrix,src_add,dst_add) |
---|
| 623 | CALL make_bicubic_remap(tabvar1(:,:,1,1),masksrc,tabinterp4d(:,:,1,1),& |
---|
| 624 | nxfin,nyfin,matrix,src_add,dst_add) |
---|
| 625 | END SELECT |
---|
| 626 | ! |
---|
| 627 | CALL Correctforconservation(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1), & |
---|
| 628 | G0%e1u,G0%e2u,G1%e1u,G1%e2u,nxfin,nyfin,posvar,imin-jpizoom+1,jmin-jpjzoom+1) |
---|
| 629 | |
---|
| 630 | ENDIF |
---|
| 631 | |
---|
| 632 | tabinterp4d(:,:,1,1) = tabinterp4d(:,:,1,1) * G1%umask(:,:,n) |
---|
| 633 | ! |
---|
| 634 | CASE('V') |
---|
| 635 | ! |
---|
| 636 | IF(MAXVAL(G1%vmask(:,:,n)) == 0) THEN |
---|
| 637 | tabinterp4d = 0.0 |
---|
| 638 | WRITE(*,*) 'only land points on level ',n |
---|
| 639 | ELSE |
---|
| 640 | ! |
---|
| 641 | |
---|
| 642 | CALL extrap_detect(G0,G1,detected_pts(:,:,n),n,'V') |
---|
| 643 | |
---|
| 644 | CALL correct_field(detected_pts(:,:,n),tabvar1,tabvar2,& |
---|
| 645 | tabvar3,G0,nav_lev,masksrc,n,'V') |
---|
| 646 | ! |
---|
| 647 | SELECT CASE(TRIM(interp_type)) |
---|
| 648 | CASE('bilinear') |
---|
| 649 | CALL get_remap_matrix(G0%gphiv,G1%gphiv, & |
---|
| 650 | G0%glamv,G1%glamv,masksrc,matrix,src_add,dst_add) |
---|
| 651 | CALL make_remap(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1),nxfin,nyfin, & |
---|
| 652 | matrix,src_add,dst_add) |
---|
| 653 | CASE('bicubic') |
---|
| 654 | CALL get_remap_bicub(G0%gphiv,G1%gphiv, & |
---|
| 655 | G0%glamv,G1%glamv,masksrc,matrix,src_add,dst_add) |
---|
| 656 | CALL make_bicubic_remap(tabvar1(:,:,1,1),masksrc,tabinterp4d(:,:,1,1),& |
---|
| 657 | nxfin,nyfin,matrix,src_add,dst_add) |
---|
| 658 | END SELECT |
---|
| 659 | ! |
---|
| 660 | CALL Correctforconservation(tabvar1(:,:,1,1),tabinterp4d(:,:,1,1), & |
---|
| 661 | G0%e1v,G0%e2v,G1%e1v,G1%e2v,nxfin,nyfin,posvar,imin-jpizoom+1,jmin-jpjzoom+1) |
---|
| 662 | |
---|
| 663 | ENDIF |
---|
| 664 | |
---|
| 665 | tabinterp4d(:,:,1,1) = tabinterp4d(:,:,1,1) * G1%vmask(:,:,n) |
---|
| 666 | ! |
---|
| 667 | END SELECT |
---|
| 668 | ! |
---|
| 669 | IF(dimg) THEN |
---|
| 670 | ! |
---|
| 671 | WRITE(inum,REC=irec) tabinterp4d(:,:,1,1) |
---|
| 672 | irec = irec + 1 |
---|
| 673 | ! |
---|
| 674 | ELSE |
---|
| 675 | ! |
---|
| 676 | dimnames(1)='x' |
---|
| 677 | dimnames(2)='y' |
---|
| 678 | IF ((iom_activated).AND.(vert_coord_name == 'z_b')) THEN |
---|
| 679 | dimnames(3)=TRIM(timedimname) |
---|
| 680 | ELSE |
---|
| 681 | dimnames(3)=vert_coord_name |
---|
| 682 | dimnames(4)=TRIM(timedimname) |
---|
| 683 | ENDIF |
---|
| 684 | ! |
---|
| 685 | IF(TRIM(varname)=='gcx' .OR. TRIM(varname)=='gcxb') THEN |
---|
| 686 | PRINT*,TRIM(varname),MAXVAL(tabinterp4d),MINVAL(tabinterp4d) |
---|
| 687 | ENDIF |
---|
| 688 | IF ((iom_activated).AND.(vert_coord_name == 'z_b')) THEN |
---|
| 689 | ALLOCATE(tabvar3d(SIZE(tabinterp4d,1),SIZE(tabinterp4d,2),SIZE(tabinterp4d,3))) |
---|
| 690 | tabvar3d=tabinterp4d(:,:,:,1) |
---|
| 691 | CALL Write_Ncdf_var(TRIM(varname),dimnames, & |
---|
| 692 | Child_file,tabvar3d,t,'double') |
---|
| 693 | DEALLOCATE(tabvar3d) |
---|
| 694 | ELSE |
---|
| 695 | CALL Write_Ncdf_var(TRIM(varname),dimnames, & |
---|
| 696 | Child_file,tabinterp4d,t,n,'double') |
---|
| 697 | ENDIF |
---|
| 698 | ! |
---|
| 699 | CALL Copy_Ncdf_att(TRIM(varname),TRIM(restart_file),Child_file) |
---|
| 700 | ! |
---|
| 701 | ENDIF |
---|
| 702 | ! |
---|
| 703 | IF(ASSOCIATED(matrix)) DEALLOCATE(matrix,src_add,dst_add) |
---|
| 704 | ! |
---|
| 705 | END DO |
---|
| 706 | ! |
---|
| 707 | DEALLOCATE(detected_pts) |
---|
| 708 | DEALLOCATE(tabinterp4d) |
---|
| 709 | DEALLOCATE(tabvar1,tabvar2,tabvar3) |
---|
| 710 | DEALLOCATE(masksrc) |
---|
| 711 | ! |
---|
| 712 | ENDIF |
---|
| 713 | |
---|
| 714 | END DO |
---|
| 715 | ! |
---|
| 716 | ! |
---|
| 717 | ! |
---|
| 718 | ! |
---|
| 719 | ! |
---|
| 720 | ! |
---|
| 721 | ! |
---|
| 722 | ! |
---|
| 723 | ! |
---|
| 724 | ALLOCATE(rotn(nxfin,nyfin,z,1),rotb(nxfin,nyfin,z,1)) |
---|
| 725 | ALLOCATE(hdivn(nxfin,nyfin,z,1),hdivb(nxfin,nyfin,z,1)) |
---|
| 726 | ! |
---|
| 727 | ! |
---|
| 728 | IF(dimg) THEN |
---|
| 729 | ALLOCATE(un(nxfin,nyfin,z,1),ub(nxfin,nyfin,z,1)) |
---|
| 730 | ALLOCATE(vn(nxfin,nyfin,z,1),vb(nxfin,nyfin,z,1)) |
---|
| 731 | irec = 6 * z + 2 |
---|
| 732 | CALL read_dimg_var(inum,irec,un,z) |
---|
| 733 | irec = 2 |
---|
| 734 | CALL read_dimg_var(inum,irec,ub,z) |
---|
| 735 | irec = 7 * z + 2 |
---|
| 736 | CALL read_dimg_var(inum,irec,vn,z) |
---|
| 737 | irec = 2 + z |
---|
| 738 | CALL read_dimg_var(inum,irec,vb,z) |
---|
| 739 | ELSE |
---|
| 740 | CALL Read_Ncdf_var('un',Child_file,un) |
---|
| 741 | CALL Read_Ncdf_var('vn',Child_file,vn) |
---|
| 742 | ! |
---|
| 743 | CALL Read_Ncdf_var('ub',Child_file,ub) |
---|
| 744 | CALL Read_Ncdf_var('vb',Child_file,vb) |
---|
| 745 | ENDIF |
---|
| 746 | ! |
---|
| 747 | |
---|
| 748 | IF(rhot == 1) THEN |
---|
| 749 | WRITE(*,*) '' |
---|
| 750 | WRITE(*,*) 'no time interpolation (time refinement ratio = 1)' |
---|
| 751 | ELSE |
---|
| 752 | WRITE(*,*) '' |
---|
| 753 | WRITE(*,*) 'time interpolation' |
---|
| 754 | now_wght = (rhot-1.)/rhot |
---|
| 755 | before_wght = 1./rhot |
---|
| 756 | ! |
---|
| 757 | ub = now_wght*un + before_wght*ub |
---|
| 758 | vb = now_wght*vn + before_wght*vb |
---|
| 759 | !----------------------- |
---|
| 760 | IF(dimg) THEN |
---|
| 761 | ALLOCATE(tn(nxfin,nyfin,z,1),tb(nxfin,nyfin,z,1)) |
---|
| 762 | irec = 8 * z + 2 |
---|
| 763 | CALL read_dimg_var(inum,irec,tn,z) |
---|
| 764 | irec = 2 * z + 2 |
---|
| 765 | CALL read_dimg_var(inum,irec,tb,z) |
---|
| 766 | ELSE |
---|
| 767 | CALL Read_Ncdf_var('tb',Child_file,tb) |
---|
| 768 | CALL Read_Ncdf_var('tn',Child_file,tn) |
---|
| 769 | ENDIF |
---|
| 770 | !---------------------- |
---|
| 771 | tb = now_wght*tn + before_wght*tb |
---|
| 772 | |
---|
| 773 | IF(dimg) THEN |
---|
| 774 | irec = 2 * z + 2 |
---|
| 775 | CALL write_dimg_var(inum,irec,tb,z) |
---|
| 776 | ELSE |
---|
| 777 | dimnames(1)='x' |
---|
| 778 | dimnames(2)='y' |
---|
| 779 | dimnames(3)='z' |
---|
| 780 | dimnames(4)=TRIM(timedimname) |
---|
| 781 | CALL Write_Ncdf_var('tb',dimnames,Child_file,tb,'double') |
---|
| 782 | ENDIF |
---|
| 783 | ! |
---|
| 784 | DEALLOCATE(tn,tb) |
---|
| 785 | !---------------------- |
---|
| 786 | IF(dimg) THEN |
---|
| 787 | ALLOCATE(sn(nxfin,nyfin,z,1),sb(nxfin,nyfin,z,1)) |
---|
| 788 | irec = 9 * z + 2 |
---|
| 789 | CALL read_dimg_var(inum,irec,sn,z) |
---|
| 790 | irec = 3 * z + 2 |
---|
| 791 | CALL read_dimg_var(inum,irec,sb,z) |
---|
| 792 | ELSE |
---|
| 793 | CALL Read_Ncdf_var('sb',Child_file,sb) |
---|
| 794 | CALL Read_Ncdf_var('sn',Child_file,sn) |
---|
| 795 | ENDIF |
---|
| 796 | !---------------------- |
---|
| 797 | sb = now_wght*sn + before_wght*sb |
---|
| 798 | |
---|
| 799 | IF(dimg) THEN |
---|
| 800 | irec = 3 * z + 2 |
---|
| 801 | CALL write_dimg_var(inum,irec,sb,z) |
---|
| 802 | ELSE |
---|
| 803 | dimnames(1)='x' |
---|
| 804 | dimnames(2)='y' |
---|
| 805 | dimnames(3)='z' |
---|
| 806 | dimnames(4)=TRIM(timedimname) |
---|
| 807 | CALL Write_Ncdf_var('sb',dimnames,Child_file,sb,'double') |
---|
| 808 | ENDIF |
---|
| 809 | !---------------------- |
---|
| 810 | DEALLOCATE(sn,sb) |
---|
| 811 | ! |
---|
| 812 | IF(dimg) THEN |
---|
| 813 | ALLOCATE(gcx(nxfin,nyfin,1,1),gcxb(nxfin,nyfin,1,1)) |
---|
| 814 | irec = 12 * z + 2 |
---|
| 815 | CALL read_dimg_var(inum,irec,gcx,1) |
---|
| 816 | irec = 12 * z + 3 |
---|
| 817 | CALL read_dimg_var(inum,irec,gcxb,1) |
---|
| 818 | ELSE |
---|
| 819 | CALL Read_Ncdf_var('gcx',Child_file,gcx) |
---|
| 820 | CALL Read_Ncdf_var('gcxb',Child_file,gcxb) |
---|
| 821 | ENDIF |
---|
| 822 | !---------------------- |
---|
| 823 | gcxb = now_wght*gcx + before_wght*gcxb |
---|
| 824 | |
---|
| 825 | IF(dimg) THEN |
---|
| 826 | irec = 12 * z + 3 |
---|
| 827 | CALL write_dimg_var(inum,irec,gcxb,1) |
---|
| 828 | ELSE |
---|
| 829 | dimnames(1)='x' |
---|
| 830 | dimnames(2)='y' |
---|
| 831 | dimnames(3)='z_b' |
---|
| 832 | dimnames(4)=TRIM(timedimname) |
---|
| 833 | CALL Write_Ncdf_var('gcxb',dimnames,Child_file,gcxb,'double') |
---|
| 834 | ENDIF |
---|
| 835 | !----------------------- |
---|
| 836 | PRINT*,' gcx = ',MAXVAL(gcx),MINVAL(gcx) |
---|
| 837 | PRINT*,' gcxb = ',MAXVAL(gcxb),MINVAL(gcxb) |
---|
| 838 | |
---|
| 839 | DEALLOCATE(gcx,gcxb) |
---|
| 840 | ! |
---|
| 841 | IF(dimg) THEN |
---|
| 842 | ALLOCATE(sshn(nxfin,nyfin,1,1),sshb(nxfin,nyfin,1,1)) |
---|
| 843 | irec = 12 * z + 5 |
---|
| 844 | CALL read_dimg_var(inum,irec,sshn,1) |
---|
| 845 | irec = 12 * z + 4 |
---|
| 846 | CALL read_dimg_var(inum,irec,sshb,1) |
---|
| 847 | ELSE |
---|
| 848 | CALL Read_Ncdf_var('sshb',Child_file,sshb) |
---|
| 849 | CALL Read_Ncdf_var('sshn',Child_file,sshn) |
---|
| 850 | ENDIF |
---|
| 851 | !---------------------- |
---|
| 852 | sshb = now_wght*sshn + before_wght*sshb |
---|
| 853 | |
---|
| 854 | IF(dimg) THEN |
---|
| 855 | irec = 12 * z + 4 |
---|
| 856 | CALL write_dimg_var(inum,irec,sshb,1) |
---|
| 857 | ELSE |
---|
| 858 | dimnames(1)='x' |
---|
| 859 | dimnames(2)='y' |
---|
| 860 | dimnames(3)='z_b' |
---|
| 861 | dimnames(4)=TRIM(timedimname) |
---|
| 862 | CALL Write_Ncdf_var('sshb',dimnames,Child_file,sshb,'double') |
---|
| 863 | ENDIF |
---|
| 864 | ! |
---|
| 865 | DEALLOCATE(sshb,sshn) |
---|
| 866 | |
---|
| 867 | ! |
---|
| 868 | ENDIF |
---|
| 869 | ! |
---|
| 870 | WRITE(*,*) 'Compute hdivn,rotn with new interpolated fields ...' |
---|
| 871 | ! |
---|
| 872 | !fmask:land/ocean mask at f-point (=0. or 1.)=shlat along lateral boundaries |
---|
| 873 | ! |
---|
| 874 | ALLOCATE(fse3u(nxfin,nyfin,z),fse3v(nxfin,nyfin,z),fse3t(nxfin,nyfin,z)) |
---|
| 875 | ! |
---|
| 876 | IF(partial_steps) THEN |
---|
| 877 | status = Read_Bathymeter(TRIM(Childbathymeter),G1) |
---|
| 878 | CALL get_scale_factors( G1,fse3t,fse3u,fse3v ) |
---|
| 879 | ELSE |
---|
| 880 | fse3t(:,:,:) = 1.0 |
---|
| 881 | fse3u(:,:,:) = 1.0 |
---|
| 882 | fse3v(:,:,:) = 1.0 |
---|
| 883 | ENDIF |
---|
| 884 | ! |
---|
| 885 | ! |
---|
| 886 | DO k = 1,z |
---|
| 887 | ! |
---|
| 888 | DO j = 2,nyfin-1 |
---|
| 889 | DO i = 2, nxfin-1 |
---|
| 890 | ! |
---|
| 891 | hdivn(i,j,k,1) = (G1%e2u(i,j)*fse3u(i,j,k)*un(i,j,k,1)-G1%e2u(i-1,j)*fse3u(i-1,j,k)*un(i-1,j,k,1) & |
---|
| 892 | + G1%e1v(i,j)*fse3v(i,j,k)*vn(i,j,k,1)-G1%e1v(i,j-1)*fse3v(i,j-1,k)*vn(i,j-1,k,1)) & |
---|
| 893 | / ( G1%e1t(i,j) * G1%e2t(i,j) * fse3t(i,j,k) ) |
---|
| 894 | ! |
---|
| 895 | hdivb(i,j,k,1) = (G1%e2u(i,j)*fse3u(i,j,k)*ub(i,j,k,1)-G1%e2u(i-1,j)*fse3u(i-1,j,k)*ub(i-1,j,k,1) & |
---|
| 896 | + G1%e1v(i,j)*fse3v(i,j,k)*vb(i,j,k,1)-G1%e1v(i,j-1)*fse3v(i,j-1,k)*vb(i,j-1,k,1)) & |
---|
| 897 | / ( G1%e1t(i,j) * G1%e2t(i,j) * fse3t(i,j,k) ) |
---|
| 898 | ! |
---|
| 899 | END DO |
---|
| 900 | END DO |
---|
| 901 | ! |
---|
| 902 | ! |
---|
| 903 | hdivn(1:2,:,k,1) = 0. |
---|
| 904 | hdivn(:,1:2,k,1) = 0. |
---|
| 905 | hdivn(nxfin-1:nxfin,:,k,1) = 0. |
---|
| 906 | hdivn(:,nyfin-1:nyfin,k,1) = 0. |
---|
| 907 | hdivb(1:2,:,k,1) = 0. |
---|
| 908 | hdivb(:,1:2,k,1) = 0. |
---|
| 909 | hdivb(nxfin-1:nxfin,:,k,1) = 0. |
---|
| 910 | hdivb(:,nyfin-1:nyfin,k,1) = 0. |
---|
| 911 | ! |
---|
| 912 | DO j = 1, nyfin-1 |
---|
| 913 | DO i = 1, nxfin-1 |
---|
| 914 | ! |
---|
| 915 | rotn(i,j,k,1) = (G1%e2v(i+1,j)*vn(i+1,j,k,1)-G1%e2v(i,j)*vn(i,j,k,1) & |
---|
| 916 | - G1%e1u(i,j+1)*un(i,j+1,k,1)+G1%e1u(i,j)*un(i,j,k,1)) & |
---|
| 917 | * G1%fmask(i,j,k) / ( G1%e1f(i,j) * G1%e2f(i,j)) |
---|
| 918 | ! |
---|
| 919 | rotb(i,j,k,1) = (G1%e2v(i+1,j)*vb(i+1,j,k,1)-G1%e2v(i,j)*vb(i,j,k,1) & |
---|
| 920 | - G1%e1u(i,j+1)*ub(i,j+1,k,1)+G1%e1u(i,j)*ub(i,j,k,1)) & |
---|
| 921 | * G1%fmask(i,j,k) / ( G1%e1f(i,j) * G1%e2f(i,j)) |
---|
| 922 | ! |
---|
| 923 | END DO |
---|
| 924 | END DO |
---|
| 925 | ! |
---|
| 926 | END DO |
---|
| 927 | ! |
---|
| 928 | PRINT*,' hdivn = ',MAXVAL(hdivn(2:nxfin-1,2:nyfin-1,:,1)),MINVAL(hdivn(2:nxfin-1,2:nyfin-1,:,1)) |
---|
| 929 | PRINT*,' hdivb = ',MAXVAL(hdivb(2:nxfin-1,2:nyfin-1,:,1)),MINVAL(hdivb(2:nxfin-1,2:nyfin-1,:,1)) |
---|
| 930 | PRINT*,' rotn = ',MAXVAL(rotn(2:nxfin-1,2:nyfin-1,:,1)),MINVAL(rotn(2:nxfin-1,2:nyfin-1,:,1)) |
---|
| 931 | PRINT*,' rotb = ',MAXVAL(rotb(2:nxfin-1,2:nyfin-1,:,1)),MINVAL(rotb(2:nxfin-1,2:nyfin-1,:,1)) |
---|
| 932 | ! |
---|
| 933 | IF(dimg) THEN |
---|
| 934 | ! |
---|
| 935 | irec = 10 * z + 2 |
---|
| 936 | DO k=1,z |
---|
| 937 | WRITE(inum,REC=irec) rotn(:,:,k,1) |
---|
| 938 | irec = irec+1 |
---|
| 939 | ENDDO |
---|
| 940 | ! |
---|
| 941 | irec = 4 * z + 2 |
---|
| 942 | DO k=1,z |
---|
| 943 | WRITE(inum,REC=irec) rotb(:,:,k,1) |
---|
| 944 | irec = irec+1 |
---|
| 945 | ENDDO |
---|
| 946 | ! |
---|
| 947 | irec = 11 * z + 2 |
---|
| 948 | DO k=1,z |
---|
| 949 | WRITE(inum,REC=irec) hdivn(:,:,k,1) |
---|
| 950 | irec = irec+1 |
---|
| 951 | ENDDO |
---|
| 952 | ! |
---|
| 953 | irec = 5 * z + 2 |
---|
| 954 | DO k=1,z |
---|
| 955 | WRITE(inum,REC=irec) hdivb(:,:,k,1) |
---|
| 956 | irec = irec+1 |
---|
| 957 | ENDDO |
---|
| 958 | ! |
---|
| 959 | CLOSE(inum) |
---|
| 960 | ELSE |
---|
| 961 | dimnames(1)='x' |
---|
| 962 | dimnames(2)='y' |
---|
| 963 | dimnames(3)='z' |
---|
| 964 | dimnames(4)=TRIM(timedimname) |
---|
| 965 | CALL Write_Ncdf_var('rotn',dimnames,Child_file,rotn,'double') |
---|
| 966 | CALL Copy_Ncdf_att('rotn',TRIM(restart_file),Child_file) |
---|
| 967 | CALL Write_Ncdf_var('hdivn',dimnames,Child_file,hdivn,'double') |
---|
| 968 | CALL Copy_Ncdf_att('hdivn',TRIM(restart_file),Child_file) |
---|
| 969 | CALL Write_Ncdf_var('rotb',dimnames,Child_file,rotb,'double') |
---|
| 970 | CALL Copy_Ncdf_att('rotb',TRIM(restart_file),Child_file) |
---|
| 971 | CALL Write_Ncdf_var('hdivb',dimnames,Child_file,hdivb,'double') |
---|
| 972 | CALL Copy_Ncdf_att('hdivb',TRIM(restart_file),Child_file) |
---|
| 973 | ENDIF |
---|
| 974 | ! |
---|
| 975 | WRITE(*,*) ' ' |
---|
| 976 | WRITE(*,*) '******* restart file successfully created *******' |
---|
| 977 | WRITE(*,*) ' ' |
---|
| 978 | ! |
---|
| 979 | STOP |
---|
| 980 | END PROGRAM create_rstrt |
---|
| 981 | |
---|
| 982 | |
---|