Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
No results found
Show changes
Showing
with 480 additions and 160 deletions
......@@ -183,6 +183,9 @@
!-----------------------------------------------------------------------
&namzdf ! vertical physics (default: NO selection)
!-----------------------------------------------------------------------
! ! adaptive-implicit vertical advection
ln_zad_Aimp = .true. ! Courant number dependent scheme (Shchepetkin 2015)
!
! ! type of vertical closure (required)
ln_zdfcst = .true. ! constant mixing
ln_zdfric = .false. ! local Richardson dependent formulation (T => fill namzdf_ric)
......
......@@ -183,6 +183,9 @@
!-----------------------------------------------------------------------
&namzdf ! vertical physics (default: NO selection)
!-----------------------------------------------------------------------
! ! adaptive-implicit vertical advection
ln_zad_Aimp = .true. ! Courant number dependent scheme (Shchepetkin 2015)
!
! ! type of vertical closure (required)
ln_zdfcst = .true. ! constant mixing
ln_zdfric = .false. ! local Richardson dependent formulation (T => fill namzdf_ric)
......
......@@ -183,6 +183,9 @@
!-----------------------------------------------------------------------
&namzdf ! vertical physics (default: NO selection)
!-----------------------------------------------------------------------
! ! adaptive-implicit vertical advection
ln_zad_Aimp = .true. ! Courant number dependent scheme (Shchepetkin 2015)
!
! ! type of vertical closure (required)
ln_zdfcst = .true. ! constant mixing
ln_zdfric = .false. ! local Richardson dependent formulation (T => fill namzdf_ric)
......
......@@ -183,6 +183,9 @@
!-----------------------------------------------------------------------
&namzdf ! vertical physics (default: NO selection)
!-----------------------------------------------------------------------
! ! adaptive-implicit vertical advection
ln_zad_Aimp = .true. ! Courant number dependent scheme (Shchepetkin 2015)
!
! ! type of vertical closure (required)
ln_zdfcst = .true. ! constant mixing
ln_zdfric = .false. ! local Richardson dependent formulation (T => fill namzdf_ric)
......
......@@ -183,6 +183,9 @@
!-----------------------------------------------------------------------
&namzdf ! vertical physics (default: NO selection)
!-----------------------------------------------------------------------
! ! adaptive-implicit vertical advection
ln_zad_Aimp = .true. ! Courant number dependent scheme (Shchepetkin 2015)
!
! ! type of vertical closure (required)
ln_zdfcst = .true. ! constant mixing
ln_zdfric = .false. ! local Richardson dependent formulation (T => fill namzdf_ric)
......
......@@ -92,16 +92,7 @@ CONTAINS
DO_2D( 0, 0, 0, 0 )
zhu(ji,jj) = 0.5_wp * ( zht(ji,jj) + zht(ji+1,jj) )
END_2D
CALL lbc_lnk( 'usrdef_zgr', zhu, 'U', 1. ) ! boundary condition: this mask the surrouding grid-points
! ! ==>>> set by hand non-zero value on first/last columns & rows
DO ji = mi0(1), mi1(1) ! first row of global domain only
zhu(ji,2) = zht(ji,2)
END DO
DO ji = mi0(jpiglo), mi1(jpiglo) ! last row of global domain only
zhu(ji,2) = zht(ji,2)
END DO
zhu(:,1) = zhu(:,2)
zhu(:,3) = zhu(:,2)
CALL lbc_lnk( 'usrdef_zgr', zhu, 'U', 1._wp ) ! boundary condition: this mask the surrouding grid-points
!
CALL zgr_z1d( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d ) ! Reference z-coordinate system
!
......@@ -111,7 +102,7 @@ CONTAINS
! no ocean cavities : top ocean level is ONE, except over land
! the ocean basin surrounded by land (1+nn_hls grid-points) set through lbc_lnk call
z2d(:,:) = 1._wp ! surface ocean is the 1st level
CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1. ) ! closed basin, see userdef_nam.F90
CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1._wp ) ! closed basin, see userdef_nam.F90
k_top(:,:) = NINT( z2d(:,:) )
!
!
......
......@@ -101,7 +101,7 @@ CONTAINS
taum(ji,jj) = zmod
wndm(ji,jj) = SQRT( zmod * zcoef )
END_2D
CALL lbc_lnk( 'usrdef_sbc', taum(:,:), 'T', 1. , wndm(:,:), 'T', 1. )
CALL lbc_lnk( 'usrdef_sbc', taum(:,:), 'T', 1._wp , wndm(:,:), 'T', 1._wp )
!
END SUBROUTINE usrdef_sbc_oce
......
......@@ -191,7 +191,7 @@ CONTAINS
!
z2d(:,:) = REAL( jpkm1 , wp ) ! flat bottom
!
CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1. ) ! set surrounding land to zero (closed boundaries)
CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1._wp ) ! set surrounding land to zero (closed boundaries)
!
zylim0 = 10000._wp ! +10km
zylim1 = 2010000._wp ! 2010km
......
......@@ -103,7 +103,7 @@ CONTAINS
END_2D
!
CALL lbc_lnk('usrdef_istate', pssh, 'T', 1. ) ! apply boundary conditions
CALL lbc_lnk('usrdef_istate', pssh, 'T', 1._wp ) ! apply boundary conditions
!
END SUBROUTINE usr_def_istate_ssh
......
......@@ -122,7 +122,7 @@ CONTAINS
END DO
END_2D
!
CALL lbc_lnk( 'usrdef_istate', pu, 'U', -1., pv, 'V', -1. )
CALL lbc_lnk( 'usrdef_istate', pu, 'U', -1._wp, pv, 'V', -1._wp )
!
END SUBROUTINE usr_def_istate
......
......@@ -189,7 +189,7 @@ CONTAINS
!
z2d(:,:) = REAL( jpkm1 , wp ) ! flat bottom
!
CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1. ) ! set surrounding land to zero (closed boundaries)
CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1._wp ) ! set surrounding land to zero (closed boundaries)
!
k_bot(:,:) = NINT( z2d(:,:) ) ! =jpkm1 over the ocean point, =0 elsewhere
!
......
......@@ -58,7 +58,6 @@
!-----------------------------------------------------------------------
&namwad ! Wetting and Drying (WaD) (default: OFF)
!-----------------------------------------------------------------------
ln_wd_il = .false ! T/F activation of iterative limiter
ln_wd_dl = .true. ! T/F activation of directional limiter
ln_wd_dl_bc = .true. ! T/F Directional limiteer Baroclinic option
ln_wd_dl_rmp = .true. ! T/F Turn on directional limiter ramp
......@@ -66,7 +65,6 @@
rn_wdmin1 = 0.2 ! Minimum wet depth on dried cells
rn_wdmin2 = 0.0001 ! Tolerance of min wet depth on dried cells
rn_wdld = 2.5 ! Land elevation below which WaD is allowed
nn_wdit = 20 ! Max iterations for WaD limiter
/
!!======================================================================
!! *** Surface Boundary Condition namelists *** !!
......@@ -156,10 +154,7 @@
cn_ice = 'none' !
nn_ice_dta = 0 ! = 0, bdy data are equal to the initial state
! = 1, bdy data are read in 'bdydata .nc' files
rn_ice_tem = 270. ! si3 only: arbitrary temperature of incoming sea ice
rn_ice_sal = 10. ! si3 only: -- salinity --
rn_ice_age = 30. ! si3 only: -- age --
!
ln_tra_dmp =.false. ! open boudaries conditions for tracers
ln_dyn3d_dmp =.false. ! open boundary condition for baroclinic velocities
rn_time_dmp = 1. ! Damping time scale in days
......
!-----------------------------------------------------------------------
&namwad ! Wetting and drying
!-----------------------------------------------------------------------
ln_wd_il = .false ! T/F activation of iterative limiter for wetting and drying scheme
ln_wd_dl = .true. ! T/F activation of directional llimiter for wetting drying scheme
ln_wd_dl_bc = .true. ! T/F Directional limiteer Baroclinic option
ln_wd_dl_rmp = .true. ! T/F Turn on directional limiter ramp
......@@ -9,7 +8,6 @@
rn_wdmin1 = 0.2 ! Minimum wet depth on dried cells
rn_wdmin2 = 0.0001 ! Tolerance of min wet depth on dried cells
rn_wdld = 2.5 ! Land elevation below which wetting/drying is allowed
nn_wdit = 20 ! Max iterations for W/D limiter
rn_wd_sbcdep = 5.0 ! Depth at which to taper sbc fluxes
rn_wd_sbcfra = 0.999 ! Fraction of SBC fluxes at taper depth (Must be <1)
/
......@@ -44,7 +44,7 @@ CONTAINS
!! ** Purpose : Initialization of the dynamics and tracers
!! Here WAD_TEST_CASES configuration
!!
q !! ** Method : - set temprature field
!! ** Method : - set temprature field
!! - set salinity field
!!----------------------------------------------------------------------
REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pdept ! depth of t-point [m]
......
......@@ -232,7 +232,7 @@ CONTAINS
DO_2D( 0, 0, 0, 0 )
zhu(ji,jj) = 0.5_wp * ( zht(ji,jj) + zht(ji+1,jj) )
END_2D
CALL lbc_lnk( 'usrdef_zgr', zhu, 'U', 1. ) ! boundary condition: this mask the surrounding grid-points
CALL lbc_lnk( 'usrdef_zgr', zhu, 'U', 1._wp ) ! boundary condition: this mask the surrounding grid-points
! ! ==>>> set by hand non-zero value on first/last columns & rows
DO ji = mi0(1), mi1(1) ! first row of global domain only
zhu(ji,:) = zht(1,:)
......@@ -245,7 +245,7 @@ CONTAINS
DO_2D( 0, 0, 0, 0 )
zhv(ji,jj) = 0.5_wp * ( zht(ji,jj) + zht(ji,jj+1) )
END_2D
CALL lbc_lnk( 'usrdef_zgr', zhv, 'V', 1. ) ! boundary condition: this mask the surrounding grid-points
CALL lbc_lnk( 'usrdef_zgr', zhv, 'V', 1._wp ) ! boundary condition: this mask the surrounding grid-points
DO jj = mj0(1), mj1(1) ! first row of global domain only
zhv(:,jj) = zht(:,jj)
END DO
......@@ -266,7 +266,7 @@ CONTAINS
z2d(:,mj0(1):mj1(1)) = 0._wp
z2d(:,mj0(jpjglo):mj1(jpjglo)) = 0._wp
CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1. ) ! closed basin, see userdef_nam.F90
CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1._wp ) ! closed basin, see userdef_nam.F90
k_top(:,:) = NINT( z2d(:,:) )
!
!
......@@ -301,15 +301,15 @@ CONTAINS
pe3vw(ji,jj,jk) = zwet * z1_jpkm1
END DO
END_2D
CALL lbc_lnk( 'usrdef_zgr', pdept, 'T', 1. )
CALL lbc_lnk( 'usrdef_zgr', pdepw, 'T', 1. )
CALL lbc_lnk( 'usrdef_zgr', pe3t , 'T', 1. )
CALL lbc_lnk( 'usrdef_zgr', pe3w , 'T', 1. )
CALL lbc_lnk( 'usrdef_zgr', pe3u , 'U', 1. )
CALL lbc_lnk( 'usrdef_zgr', pe3uw, 'U', 1. )
CALL lbc_lnk( 'usrdef_zgr', pe3f , 'F', 1. )
CALL lbc_lnk( 'usrdef_zgr', pe3v , 'V', 1. )
CALL lbc_lnk( 'usrdef_zgr', pe3vw, 'V', 1. )
CALL lbc_lnk( 'usrdef_zgr', pdept, 'T', 1._wp )
CALL lbc_lnk( 'usrdef_zgr', pdepw, 'T', 1._wp )
CALL lbc_lnk( 'usrdef_zgr', pe3t , 'T', 1._wp )
CALL lbc_lnk( 'usrdef_zgr', pe3w , 'T', 1._wp )
CALL lbc_lnk( 'usrdef_zgr', pe3u , 'U', 1._wp )
CALL lbc_lnk( 'usrdef_zgr', pe3uw, 'U', 1._wp )
CALL lbc_lnk( 'usrdef_zgr', pe3f , 'F', 1._wp )
CALL lbc_lnk( 'usrdef_zgr', pe3v , 'V', 1._wp )
CALL lbc_lnk( 'usrdef_zgr', pe3vw, 'V', 1._wp )
WHERE( pe3t (:,:,:) == 0._wp ) pe3t (:,:,:) = 1._wp
WHERE( pe3u (:,:,:) == 0._wp ) pe3u (:,:,:) = 1._wp
WHERE( pe3v (:,:,:) == 0._wp ) pe3v (:,:,:) = 1._wp
......
......@@ -748,7 +748,7 @@
EXTERNAL :: init_e2t, init_e2u, init_e2v, init_e2f
ln_perio=.FALSE.
if( jperio ==1 .OR. jperio==2 .OR. jperio==4) ln_perio=.TRUE.
IF( jperio == 1 .OR. jperio == 2 .OR. jperio == 4 .OR. jperio == 6 ) ln_perio=.TRUE.
CALL Agrif_Init_variable(e1t_id, procname = init_e1t)
CALL Agrif_Init_variable(e1u_id, procname = init_e1u)
......@@ -1216,8 +1216,8 @@
CALL ctl_stop( 'STOP', 'AGRIF zoom imin must be < imax' )
ENDIF
IF ( (Agrif_Parent(jperio)==4).OR.(Agrif_Parent(jperio)==1) ) THEN
IF ( (jperio==4).OR.(jperio==1) ) THEN ! Cyclic east-west zoom
IF ( Agrif_Parent(jperio) == 1 .OR. Agrif_Parent(jperio) == 4 .OR. Agrif_Parent(jperio) == 6 ) THEN
IF ( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 ) THEN ! Cyclic east-west zoom
lk_west = .FALSE. ; lk_east = .FALSE.
! Checks:
IF ( imin/=1-Agrif_Parent(nbghostcells_x_w) ) THEN
......@@ -1258,8 +1258,8 @@
CALL ctl_stop( 'STOP', 'AGRIF zoom jmin must be < jmax' )
ENDIF
IF ( Agrif_Parent(jperio)==4 ) THEN
IF (jperio==4) THEN ! North-Fold
IF ( Agrif_Parent(jperio) == 4 .OR. Agrif_Parent(jperio) == 6 ) THEN
IF ( jperio == 4 .OR. jperio == 6 ) THEN ! North-Fold
lk_north = .FALSE.
! Checks:
IF ( jmax/=Agrif_Parent(Nj0glo)+1-Agrif_Parent(nbghostcells_y_s)) THEN
......@@ -1300,10 +1300,10 @@
IF (.NOT.lk_south) nbghostcells_y_s = 1
IF (.NOT.lk_north) nbghostcells_y_n = 1
IF ((jperio == 1).OR.(jperio == 4)) THEN
IF ( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 ) THEN
nbghostcells_x_w = 0 ; nbghostcells_x_e = 0
ENDIF
IF (jperio == 4) THEN
IF ( jperio == 4 .OR. jperio == 6 ) THEN
nbghostcells_y_n = 0
ENDIF
......
#!/bin/bash
#
# this sed makes sure that we alwas specify the _wp in the calls to lbc_lnk
#
for ff in $( find src -type f ) */*/MY_SRC/* # for all nemo files
do
# look for something like ", 'X', 1. )" and replace it by ", 'X', 1.wp )"
# - X can be TUVFW
# - the "1" can be followed by a "." or not. It must be followed by "," or ")" or "&"
# - any number of blanks is acceptable.
sed -e "s/\(, *'[TUVFW]' *, *-*1\)\.*\( *[,)&]\)/\1._wp\2/g" $ff > tmp
mv tmp $ff
done
......@@ -123,18 +123,20 @@ PROGRAM create_bathy
! check grids: if identical then do not interpolate
identical_grids = .FALSE.
IF( SIZE(G0%nav_lat,1) == SIZE(G1%nav_lat,1) .AND. SIZE(G0%nav_lat,2) == SIZE(G1%nav_lat,2) .AND. &
& SIZE(G0%nav_lon,1) == SIZE(G1%nav_lon,1) .AND. SIZE(G0%nav_lon,2) == SIZE(G1%nav_lon,2) ) THEN
IF( MAXVAL( ABS(G0%nav_lat(:,:)- G1%nav_lat(:,:)) ) < 0.0001 .AND. &
& MAXVAL( ABS(G0%nav_lon(:,:)- G1%nav_lon(:,:)) ) < 0.0001 ) THEN
WRITE(*,*) ''
WRITE(*,*) 'same grid between parent and child domains => NO INTERPOLATION'
WRITE(*,*) ''
G1%bathy_meter = G0%bathy_meter
identical_grids = .TRUE.
ENDIF
ENDIF
!clem: comment these lines if one wants to create a zoom 1:1 and to change the bathymetry in the zoom
!
! IF( SIZE(G0%nav_lat,1) == SIZE(G1%nav_lat,1) .AND. SIZE(G0%nav_lat,2) == SIZE(G1%nav_lat,2) .AND. &
! & SIZE(G0%nav_lon,1) == SIZE(G1%nav_lon,1) .AND. SIZE(G0%nav_lon,2) == SIZE(G1%nav_lon,2) ) THEN
! IF( MAXVAL( ABS(G0%nav_lat(:,:)- G1%nav_lat(:,:)) ) < 0.0001 .AND. &
! & MAXVAL( ABS(G0%nav_lon(:,:)- G1%nav_lon(:,:)) ) < 0.0001 ) THEN
! WRITE(*,*) ''
! WRITE(*,*) 'same grid between parent and child domains => NO INTERPOLATION'
! WRITE(*,*) ''
! G1%bathy_meter = G0%bathy_meter
! identical_grids = .TRUE.
! ENDIF
! ENDIF
IF( .NOT.new_topo ) type_bathy_interp = 2 ! only one which works
!
......@@ -270,8 +272,8 @@ PROGRAM create_bathy
ALLOCATE(masksrc(SIZE(G0%bathy_meter,1),SIZE(G0%bathy_meter,2)))
ALLOCATE(bathy_test(nxfin,nyfin))
!
WHERE(G0%bathy_meter.LE.0) ; masksrc = .FALSE. ;
ELSEWHERE ; masksrc = .TRUE. ;
WHERE(G0%bathy_meter.LE.0.) ; masksrc = .FALSE. ;
ELSEWHERE ; masksrc = .TRUE. ;
END WHERE
!
! compute remapping matrix thanks to SCRIP package
......
......@@ -99,7 +99,7 @@ PROGRAM create_rstrt
! create this file
!
CALL set_child_name(restart_file,Child_file)
status = nf90_create(Child_file,NF90_WRITE,ncid)
status = nf90_create(Child_file,NF90_NETCDF4,ncid)
status = nf90_close(ncid)
WRITE(*,*) 'Child grid restart file name = ',TRIM(Child_file)
WRITE(*,*) ''
......
......@@ -210,21 +210,21 @@ CONTAINS
SELECT CASE(typevar)
CASE('T')
IF(MOD(irafx,2)==1) THEN ! odd
zx = 1 ; zy = 1 ; jdecx = FLOOR(irafx/2.) ; jdecy = FLOOR(irafy/2.)
zx = 1 ; zy = 1 ; jdecx = FLOOR(REAL(irafx)/2.) ; jdecy = FLOOR(REAL(irafy)/2.)
ELSE ! even
zx = 2 ; zy = 2 ; jdecx = FLOOR(irafx/2.) ; jdecy = FLOOR(irafy/2.)
zx = 2 ; zy = 2 ; jdecx = FLOOR(REAL(irafx)/2.) ; jdecy = FLOOR(REAL(irafy)/2.)
ENDIF
CASE('U')
IF(MOD(irafx,2)==1) THEN ! odd
zx = 1 ; zy = 1 ; jdecx = irafx - 1 ; jdecy = FLOOR(irafy/2.)
zx = 1 ; zy = 1 ; jdecx = irafx - 1 ; jdecy = FLOOR(REAL(irafy)/2.)
ELSE ! even
zx = 1 ; zy = 2 ; jdecx = irafx - 1 ; jdecy = FLOOR(irafy/2.)
zx = 1 ; zy = 2 ; jdecx = irafx - 1 ; jdecy = FLOOR(REAL(irafy)/2.)
ENDIF
CASE('V')
IF(MOD(irafx,2)==1) THEN ! odd
zx = 1 ; zy = 1 ; jdecx = FLOOR(irafx/2.) ; jdecy = irafy - 1
zx = 1 ; zy = 1 ; jdecx = FLOOR(REAL(irafx)/2.) ; jdecy = irafy - 1
ELSE ! even
zx = 2 ; zy = 1 ; jdecx = FLOOR(irafx/2.) ; jdecy = irafy - 1
zx = 2 ; zy = 1 ; jdecx = FLOOR(REAL(irafx)/2.) ; jdecy = irafy - 1
ENDIF
CASE('F')
IF(MOD(irafx,2)==1) THEN ! odd
......@@ -238,21 +238,24 @@ CONTAINS
DO jj = 1, nyf
jjf = jj - jdecy
jjc = j_min + FLOOR((jjf-1.) / irafy)
jjc = j_min + FLOOR(REAL(jjf-1) / REAL(irafy))
DO ji = 1, nxf
jif = ji - jdecx
jic = i_min + FLOOR((jif-1.) / irafx)
Bx = MOD( zx*jif-1, zx*irafx ) / REAL(zx*irafx)
By = MOD( zy*jjf-1, zy*irafy ) / REAL(zy*irafy)
jic = i_min + FLOOR(REAL(jif-1) / REAL(irafx))
Bx = REAL(MODULO( zx*jif-1, zx*irafx )) / REAL(zx*irafx)
By = REAL(MODULO( zy*jjf-1, zy*irafy )) / REAL(zy*irafy)
Ax = 1. - Bx
Ay = 1. - By
jic1 = MIN( nxc, jic+1 ) ! avoid out of bounds for tabin below
jjc1 = MIN( nyc, jjc+1 ) ! --
jic = MAX( 1, jic ) ! avoid out of bounds
jjc = MAX( 1, jjc ) ! avoid out of bounds
tabout(ji,jj) = ( Bx * tabin(jic1,jjc ) + Ax * tabin(jic,jjc ) ) * Ay + &
& ( Bx * tabin(jic1,jjc1) + Ax * tabin(jic,jjc1) ) * By
END DO
......@@ -503,7 +506,7 @@ CONTAINS
REAL*8, POINTER, DIMENSION(:,:) :: lonParent,latParent => NULL()
REAL*8, POINTER, DIMENSION(:,:) :: lonChild,latChild,latlon_temp => NULL()
REAL*8, POINTER, DIMENSION(:,:,:,:) :: tabinterp4d,tabvar1,tabvar2,tabvar3 => NULL()
REAL*8, POINTER, DIMENSION(:,:,:) :: tabinterp3d,tabvar3d => NULL()
REAL*8, POINTER, DIMENSION(:,:,:) :: tabinterp3d,tabvar3d,tabvar01,tabvar02,tabvar03 => NULL()
REAL*8, POINTER, DIMENSION(:) :: timedepth_temp,depth => NULL()
REAL*8,DIMENSION(:,:),POINTER :: matrix => NULL()
INTEGER,DIMENSION(:),POINTER :: src_add,dst_add => NULL()
......@@ -515,11 +518,17 @@ CONTAINS
LOGICAL,DIMENSION(:,:),POINTER :: masksrc => NULL()
LOGICAL :: Interpolation,conservation,Pacifique,Extrapolation,land_level
!
INTEGER :: deptht,time,i,status,ncid,t,ii,j,nb,numlon,numlat
INTEGER :: deptht,time,i,status,ncid,t,ii,j,nb,numlon,numlat,unlimdimid
!
TYPE(Coordinates) :: G0,G1
!
status = nf90_open(TRIM(filename),NF90_NOWRITE,ncid)
IF (status/=nf90_noerr) THEN
WRITE(*,*)"unable to open netcdf file : ",TRIM(filename)
STOP
ENDIF
!
!*****************
!If coarse grid is masked possibility to activate an extrapolation process
!
......@@ -541,7 +550,6 @@ CONTAINS
!
CALL Read_Ncdf_dim('x',filename,numlon)
CALL Read_Ncdf_dim('y',filename,numlat)
CALL Read_Ncdf_dim('time_counter',filename,time)
IF ( Dims_Existence( 'deptht' , filename ) ) THEN
CALL Read_Ncdf_dim('deptht',filename,deptht)
ELSE IF ( Dims_Existence( 'depthu' , filename ) ) THEN
......@@ -555,7 +563,28 @@ CONTAINS
ELSE
deptht = N
ENDIF
IF ( Dims_Existence( 'time_counter' , filename ) ) THEN
CALL Read_Ncdf_dim('time_counter',filename,time)
! check that time_counter is the unlimited dim
status = nf90_inquire(ncid, unlimiteddimid = unlimdimid)
IF ( unlimdimid == -1 ) THEN
WRITE(*,*)"time_counter should be the unlimited dimension in this file : ",filename
WRITE(*,*)" use nco command: ncks -O --mk_rec_dmn time_counter ",filename," ",filename
STOP
ENDIF
ELSEIF( Dims_Existence( 'time' , filename ) ) THEN
CALL Read_Ncdf_dim('time',filename,time)
! check that time is the unlimited dim
status = nf90_inquire(ncid, unlimiteddimid = unlimdimid)
IF ( unlimdimid == -1 ) THEN
WRITE(*,*)"time should be the unlimited dimension in this file : ",filename
WRITE(*,*)" use nco command: ncks -O --mk_rec_dmn time ",filename," ",filename
STOP
ENDIF
ELSE
time = 0
ENDIF
!
! retrieve netcdf variable name
!
......@@ -579,7 +608,7 @@ CONTAINS
status = Read_Coordinates(Childcoordinates,G1,Pacifique)
!
! check consistency of informations read in namelist
!
!
IF( imax > SIZE(G0%glamt,1) .OR. jmax > SIZE(G0%glamt,2) .OR. &
imax <= imin .OR. jmax <= jmin ) THEN
WRITE(*,*) 'ERROR ***** bad child grid definition ...'
......@@ -597,7 +626,6 @@ CONTAINS
STOP
ENDIF
!
!
! Initialization of T-mask thanks to bathymetry
!
......@@ -610,7 +638,6 @@ CONTAINS
CALL Init_mask(parent_bathy_level,G0,1,1)
!
ENDIF
!
! select coordinates to use according to variable position
!
......@@ -619,23 +646,67 @@ CONTAINS
SELECT CASE(posvar)
CASE('T')
lonParent = G0%glamt
latParent = G0%gphit
lonChild = G1%glamt
latChild = G1%gphit
mask = G1%tmask
IF ( Vars_Existence( 'glamt' , parent_coordinate_file ) ) THEN
lonParent = G0%glamt
latParent = G0%gphit
ELSEIF( Vars_Existence( 'nav_lon' , parent_coordinate_file ) ) THEN
lonParent = G0%nav_lon
latParent = G0%nav_lat
ENDIF
IF ( Vars_Existence( 'glamt' , Childcoordinates ) ) THEN
lonChild = G1%glamt
latChild = G1%gphit
ELSEIF( Vars_Existence( 'nav_lon' , Childcoordinates ) ) THEN
lonChild = G1%nav_lon
latChild = G1%nav_lat
ENDIF
!!IF ( Vars_Existence( 'tmask' , Childcoordinates ) ) THEN
mask = G1%tmask
!!ELSEIF( Vars_Existence( 'top_level' , Childcoordinates ) ) THEN
!! mask = MAX( 1., G1%top_level )
!!ELSE
!! mask = 1.
!!ENDIF
CASE('U')
lonParent = G0%glamu
latParent = G0%gphiu
lonChild = G1%glamu
latChild = G1%gphiu
mask = G1%umask
IF ( Vars_Existence( 'glamu' , parent_coordinate_file ) ) THEN
lonParent = G0%glamu
latParent = G0%gphiu
ELSE
WRITE(*,*) 'ERROR ***** missing glamu/gphiu coordinates in ', TRIM(parent_coordinate_file)
STOP
ENDIF
IF ( Vars_Existence( 'glamu' , Childcoordinates ) ) THEN
lonChild = G1%glamu
latChild = G1%gphiu
ELSE
WRITE(*,*) 'ERROR ***** missing glamu/gphiu coordinates in ', TRIM(Childcoordinates)
STOP
ENDIF
!!IF ( Vars_Existence( 'umask' , Childcoordinates ) ) THEN
mask = G1%umask
!!ELSE
!! mask = 1.
!!ENDIF
CASE('V')
lonParent = G0%glamv
latParent = G0%gphiv
lonChild = G1%glamv
latChild = G1%gphiv
mask = G1%vmask
IF ( Vars_Existence( 'glamv' , parent_coordinate_file ) ) THEN
lonParent = G0%glamv
latParent = G0%gphiv
ELSE
WRITE(*,*) 'ERROR ***** missing glamv/gphiv coordinates in ', TRIM(parent_coordinate_file)
STOP
ENDIF
IF ( Vars_Existence( 'glamv' , Childcoordinates ) ) THEN
lonChild = G1%glamv
latChild = G1%gphiv
ELSE
WRITE(*,*) 'ERROR ***** missing glamv/gphiv coordinates in ', TRIM(Childcoordinates)
STOP
ENDIF
!!IF ( Vars_Existence( 'vmask' , Childcoordinates ) ) THEN
mask = G1%vmask
!!ELSE
!! mask = 1.
!!ENDIF
END SELECT
DEALLOCATE(G0%glamu,G0%glamv,G0%gphiu,G0%gphiv)
......@@ -665,8 +736,9 @@ CONTAINS
IF ( Dims_Existence( 'depthu' , filename ) ) CALL Write_Ncdf_dim('depthu',Child_file,deptht)
IF ( Dims_Existence( 'depthv' , filename ) ) CALL Write_Ncdf_dim('depthv',Child_file,deptht)
IF ( Dims_Existence( 'depthw' , filename ) ) CALL Write_Ncdf_dim('depthw',Child_file,deptht)
IF ( Dims_Existence( 'z' , filename ) ) CALL Write_Ncdf_dim('z',Child_file,deptht)
CALL Write_Ncdf_dim('time_counter',Child_file,0)
IF ( Dims_Existence( 'z' , filename ) ) CALL Write_Ncdf_dim('z' ,Child_file,deptht)
IF ( Dims_Existence( 'time_counter' , filename ) ) CALL Write_Ncdf_dim('time_counter',Child_file,time)
IF ( Dims_Existence( 'time' , filename ) ) CALL Write_Ncdf_dim('time' ,Child_file,time)
IF( deptht .NE. 1 .AND. deptht .NE. N ) THEN
WRITE(*,*) '***'
......@@ -733,6 +805,17 @@ CONTAINS
varname = TRIM(Ncdf_varname(i))
Interpolation = .FALSE.
!
!copy time from input file to output file
CASE('time')
ALLOCATE(timedepth_temp(time))
CALL Read_Ncdf_var('time',filename,timedepth_temp)
CALL Write_Ncdf_var('time','time', &
Child_file,timedepth_temp,'float')
CALL Copy_Ncdf_att('time',TRIM(filename),Child_file)
DEALLOCATE(timedepth_temp)
varname = TRIM(Ncdf_varname(i))
Interpolation = .FALSE.
!
!copy deptht from input file to output file
CASE('deptht')
ALLOCATE(depth(deptht))
......@@ -806,84 +889,305 @@ CONTAINS
!
END SELECT
! //////////////// INTERPOLATION FOR 3D VARIABLES /////////////////////////////////////
! //////////////// INTERPOLATION FOR 2D VARIABLES /////////////////////////////////////
!
IF( Interpolation .AND. Get_NbDims(TRIM(varname),TRIM(filename)) == 3 ) THEN
IF( Interpolation .AND. Get_NbDims(TRIM(varname),TRIM(filename)) == 2 ) THEN
!
ALLOCATE(detected_pts(numlon,numlat,N))
ALLOCATE(masksrc(numlon,numlat))
!
! ******************************LOOP ON TIME*******************************************
!loop on time
DO t = 1,time
!
IF(extrapolation) THEN
WRITE(*,*) 'interpolation/extrapolation ',TRIM(varname),' for time t = ',t
ELSE
WRITE(*,*) 'interpolation ',TRIM(varname),' for time t = ',t
ENDIF
!
ALLOCATE(tabvar3d(numlon,numlat,1))
ALLOCATE(tabinterp3d(nxfin,nyfin,1))
!
CALL Read_Ncdf_var(varname,filename,tabvar3d,t)
IF(extrapolation) THEN
WRITE(*,*) 'interpolation/extrapolation ',TRIM(varname)
ELSE
WRITE(*,*) 'interpolation ',TRIM(varname)
ENDIF
!
ALLOCATE(tabvar3d(numlon,numlat,1))
ALLOCATE(tabinterp3d(nxfin,nyfin,1))
!
CALL Read_Ncdf_var(varname,filename,tabvar3d)
!
! search points where extrapolation is required
!
IF(Extrapolation) THEN
WHERE( tabvar3d .GE. 1.e+20 ) tabvar3d = 0.
CALL extrap_detect(G0,G1,detected_pts(:,:,1),1)
CALL correct_field_2d(detected_pts(:,:,1),tabvar3d,G0,masksrc,'posvar')
ELSE
masksrc = .TRUE.
ENDIF
SELECT CASE(TRIM(interp_type))
CASE('bilinear')
CALL get_remap_matrix(latParent,latChild, &
lonParent,lonChild, &
masksrc,matrix,src_add,dst_add)
CASE('bicubic')
CALL get_remap_bicub(latParent,latChild, &
lonParent,lonChild, &
masksrc,matrix,src_add,dst_add)
END SELECT
!
!
SELECT CASE(TRIM(interp_type))
CASE('bilinear')
CALL make_remap(tabvar3d(:,:,1),tabinterp3d(:,:,1),nxfin,nyfin, &
matrix,src_add,dst_add)
CASE('bicubic')
CALL make_bicubic_remap(tabvar3d(:,:,1),masksrc,tabinterp3d(:,:,1),nxfin,nyfin, &
matrix,src_add,dst_add)
END SELECT
!
IF( conservation ) CALL Correctforconservation(tabvar3d(:,:,1),tabinterp3d(:,:,1), &
G0%e1t,G0%e2t,G1%e1t,G1%e2t,nxfin,nyfin,posvar,imin,jmin)
!
IF(Extrapolation) tabinterp3d(:,:,1) = tabinterp3d(:,:,1) * mask(:,:,1)
!! tabinterp3d(:,:,1) = tabinterp3d(:,:,1) * mask(:,:,1)
!
dimnames(1)='x'
dimnames(2)='y'
!
CALL Write_Ncdf_var(TRIM(varname),dimnames(1:2),&
Child_file,tabinterp3d(:,:,1),'float')
!
DEALLOCATE(tabinterp3d)
DEALLOCATE(tabvar3d)
!
DEALLOCATE(detected_pts)
IF(ASSOCIATED(matrix)) DEALLOCATE(matrix,dst_add,src_add)
DEALLOCATE( masksrc)
CALL Copy_Ncdf_att(TRIM(varname),TRIM(filename),Child_file)
!
! //////////////// INTERPOLATION FOR 3D VARIABLES /////////////////////////////////////
!
ELSEIF( Interpolation .AND. Get_NbDims(TRIM(varname),TRIM(filename)) == 3 ) THEN
!
IF( DimUnlimited_Var(TRIM(varname),TRIM(filename)) ) THEN
ALLOCATE(detected_pts(numlon,numlat,N))
ALLOCATE(masksrc(numlon,numlat))
!
! search points where extrapolation is required
!
IF(Extrapolation) THEN
WHERE( tabvar3d .GE. 1.e+20 ) tabvar3d = 0.
IF (t .EQ. 1. ) CALL extrap_detect(G0,G1,detected_pts(:,:,1),1)
CALL correct_field_2d(detected_pts(:,:,1),tabvar3d,G0,masksrc,'posvar')
ELSE
masksrc = .TRUE.
ENDIF
! ******************************LOOP ON TIME*******************************************
!loop on time
DO t = 1,time
!
IF(extrapolation) THEN
WRITE(*,*) 'interpolation/extrapolation ',TRIM(varname),' for time t = ',t
ELSE
WRITE(*,*) 'interpolation ',TRIM(varname),' for time t = ',t
ENDIF
!
ALLOCATE(tabvar3d(numlon,numlat,1))
ALLOCATE(tabinterp3d(nxfin,nyfin,1))
!
CALL Read_Ncdf_var(varname,filename,tabvar3d,t)
!
! search points where extrapolation is required
!
IF(Extrapolation) THEN
WHERE( tabvar3d .GE. 1.e+20 ) tabvar3d = 0.
IF (t .EQ. 1. ) CALL extrap_detect(G0,G1,detected_pts(:,:,1),1)
CALL correct_field_2d(detected_pts(:,:,1),tabvar3d,G0,masksrc,'posvar')
ELSE
masksrc = .TRUE.
ENDIF
IF (t.EQ.1 ) THEN
IF (t.EQ.1 ) THEN
SELECT CASE(TRIM(interp_type))
CASE('bilinear')
CALL get_remap_matrix(latParent,latChild, &
lonParent,lonChild, &
masksrc,matrix,src_add,dst_add)
CASE('bicubic')
CALL get_remap_bicub(latParent,latChild, &
lonParent,lonChild, &
masksrc,matrix,src_add,dst_add)
END SELECT
!
ENDIF
!
SELECT CASE(TRIM(interp_type))
CASE('bilinear')
CALL get_remap_matrix(latParent,latChild, &
lonParent,lonChild, &
masksrc,matrix,src_add,dst_add)
CASE('bilinear')
CALL make_remap(tabvar3d(:,:,1),tabinterp3d(:,:,1),nxfin,nyfin, &
matrix,src_add,dst_add)
CASE('bicubic')
CALL make_bicubic_remap(tabvar3d(:,:,1),masksrc,tabinterp3d(:,:,1),nxfin,nyfin, &
matrix,src_add,dst_add)
END SELECT
!
IF( conservation ) CALL Correctforconservation(tabvar3d(:,:,1),tabinterp3d(:,:,1), &
G0%e1t,G0%e2t,G1%e1t,G1%e2t,nxfin,nyfin,posvar,imin,jmin)
!
IF(Extrapolation) tabinterp3d(:,:,1) = tabinterp3d(:,:,1) * mask(:,:,1)
!! tabinterp3d(:,:,1) = tabinterp3d(:,:,1) * mask(:,:,1)
!
dimnames(1)='x'
dimnames(2)='y'
dimnames(3)='time_counter'
!
CALL Write_Ncdf_var(TRIM(varname),dimnames(1:3),&
Child_file,tabinterp3d,t,'float')
!
DEALLOCATE(tabinterp3d)
DEALLOCATE(tabvar3d)
!end loop on time
END DO
!
DEALLOCATE(detected_pts)
IF(ASSOCIATED(matrix)) DEALLOCATE(matrix,dst_add,src_add)
DEALLOCATE( masksrc)
CASE('bicubic')
CALL get_remap_bicub(latParent,latChild, &
lonParent,lonChild, &
masksrc,matrix,src_add,dst_add)
ELSE
END SELECT
!
ENDIF
!
SELECT CASE(TRIM(interp_type))
CASE('bilinear')
CALL make_remap(tabvar3d(:,:,1),tabinterp3d(:,:,1),nxfin,nyfin, &
matrix,src_add,dst_add)
CASE('bicubic')
CALL make_bicubic_remap(tabvar3d(:,:,1),masksrc,tabinterp3d(:,:,1),nxfin,nyfin, &
matrix,src_add,dst_add)
END SELECT
!
IF( conservation ) CALL Correctforconservation(tabvar3d(:,:,1),tabinterp3d(:,:,1), &
G0%e1t,G0%e2t,G1%e1t,G1%e2t,nxfin,nyfin,posvar,imin,jmin)
!
IF(Extrapolation) tabinterp3d(:,:,1) = tabinterp3d(:,:,1) * mask(:,:,1)
!
dimnames(1)='x'
dimnames(2)='y'
dimnames(3)='time_counter'
!
CALL Write_Ncdf_var(TRIM(varname),dimnames(1:3),&
Child_file,tabinterp3d,t,'float')
!
DEALLOCATE(tabinterp3d)
DEALLOCATE(tabvar3d)
!end loop on time
END DO
!
DEALLOCATE(detected_pts)
IF(ASSOCIATED(matrix)) DEALLOCATE(matrix,dst_add,src_add)
DEALLOCATE( masksrc)
IF ( Dims_Existence( 'deptht' , filename ) ) dimnames(3)='deptht'
IF ( Dims_Existence( 'depthu' , filename ) ) dimnames(3)='depthu'
IF ( Dims_Existence( 'depthv' , filename ) ) dimnames(3)='depthv'
IF ( Dims_Existence( 'depthw' , filename ) ) dimnames(3)='depthw'
IF ( Dims_Existence( 'z' , filename ) ) dimnames(3)='z'
!
! loop on vertical levels
!
DO nb = 1,deptht
!
ALLOCATE(masksrc(numlon,numlat))
ALLOCATE(detected_pts(numlon,numlat,N))
!
! point detection et level n
!
land_level = .FALSE.
IF( Extrapolation ) THEN
IF(MAXVAL(mask(:,:,nb))==0.) land_level = .TRUE.
ENDIF
IF ( land_level ) THEN
!
WRITE(*,*) 'only land points on level ',nb
ALLOCATE(tabinterp3d(nxfin,nyfin,1))
tabinterp3d = 0.e0
!
CALL Write_Ncdf_var(TRIM(varname),dimnames, &
Child_file,tabinterp3d,nb,'float')
DEALLOCATE(tabinterp3d)
!
ELSE
!
ALLOCATE(tabvar01(numlon,numlat,1)) ! level k
IF(Extrapolation) ALLOCATE(tabvar02(numlon,numlat,1)) ! level k-1
IF(Extrapolation) ALLOCATE(tabvar03(numlon,numlat,1)) ! level k-2
ALLOCATE(tabinterp3d(nxfin,nyfin,1))
!
IF(Extrapolation) THEN
IF(nb==1) THEN
CALL Read_Ncdf_var(varname,filename,tabvar01,nb)
WHERE( tabvar01 .GE. 1.e+20 ) tabvar01 = 0.
ELSE IF (nb==2) THEN
CALL Read_Ncdf_var(varname,filename,tabvar02,nb-1)
CALL Read_Ncdf_var(varname,filename,tabvar01,nb)
WHERE( tabvar01 .GE. 1.e+20 ) tabvar01 = 0.
WHERE( tabvar02 .GE. 1.e+20 ) tabvar02 = 0.
ELSE
CALL Read_Ncdf_var(varname,filename,tabvar03,nb-2)
CALL Read_Ncdf_var(varname,filename,tabvar02,nb-1)
CALL Read_Ncdf_var(varname,filename,tabvar01,nb)
WHERE( tabvar01 .GE. 1.e+20 ) tabvar01 = 0.
WHERE( tabvar02 .GE. 1.e+20 ) tabvar02 = 0.
WHERE( tabvar03 .GE. 1.e+20 ) tabvar03 = 0.
ENDIF
!
CALL extrap_detect(G0,G1,detected_pts(:,:,nb),nb)
ALLOCATE(tabvar1(numlon,numlat,1,1)) ! level k
ALLOCATE(tabvar2(numlon,numlat,1,1)) ! level k-1
ALLOCATE(tabvar3(numlon,numlat,1,1)) ! level k-2
tabvar1(:,:,1,1) = tabvar01(:,:,1)
tabvar2(:,:,1,1) = tabvar02(:,:,1)
tabvar3(:,:,1,1) = tabvar03(:,:,1)
CALL correct_field(detected_pts(:,:,nb),tabvar1,tabvar2,&
tabvar3,G0,depth,masksrc,nb,'posvar')
tabvar01(:,:,1) = tabvar1(:,:,1,1)
tabvar02(:,:,1) = tabvar2(:,:,1,1)
tabvar03(:,:,1) = tabvar3(:,:,1,1)
DEALLOCATE(tabvar1,tabvar2,tabvar3)
DEALLOCATE(tabvar02,tabvar03)
ELSE
CALL Read_Ncdf_var(varname,filename,tabvar01,nb)
IF(MAXVAL(tabvar01(:,:,1))==0.) land_level = .TRUE.
masksrc = .TRUE.
ENDIF
IF( Extrapolation ) THEN
WRITE(*,*) 'interpolation/extrapolation ',TRIM(varname),' for vertical level = ',nb
ELSE
WRITE(*,*) 'interpolation ',TRIM(varname),' for vertical level = ',nb
ENDIF
!
SELECT CASE(TRIM(interp_type))
CASE('bilinear')
CALL get_remap_matrix(latParent,latChild, &
lonParent,lonChild, &
masksrc,matrix,src_add,dst_add)
CASE('bicubic')
CALL get_remap_bicub(latParent,latChild, &
lonParent,lonChild, &
masksrc,matrix,src_add,dst_add)
!
END SELECT
!
!
SELECT CASE(TRIM(interp_type))
!
CASE('bilinear')
CALL make_remap(tabvar01(:,:,1),tabinterp3d(:,:,1),nxfin,nyfin, &
matrix,src_add,dst_add)
CASE('bicubic')
CALL make_bicubic_remap(tabvar01(:,:,1),masksrc,tabinterp3d(:,:,1),nxfin,nyfin, &
matrix,src_add,dst_add)
END SELECT
!
IF( conservation ) CALL Correctforconservation(tabvar01(:,:,1),tabinterp3d(:,:,1), &
G0%e1t,G0%e2t,G1%e1t,G1%e2t,nxfin,nyfin,posvar,imin,jmin)
!
ALLOCATE(tabinterp4d(nxfin,nyfin,1,1))
tabinterp4d(:,:,1,1) = tabinterp3d(:,:,1)
IF(Extrapolation) CALL check_extrap(G1,tabinterp4d,nb)
tabinterp3d(:,:,1) = tabinterp4d(:,:,1,1)
DEALLOCATE(tabinterp4d)
!
IF(Extrapolation) tabinterp3d(:,:,1) = tabinterp3d(:,:,1) * mask(:,:,nb)
!! tabinterp3d(:,:,1) = tabinterp3d(:,:,1) * mask(:,:,nb)
!
CALL Write_Ncdf_var(TRIM(varname),dimnames, &
Child_file,tabinterp3d,nb,'float')
!
DEALLOCATE(tabinterp3d)
DEALLOCATE(tabvar01)
!
!
ENDIF
!
IF(ASSOCIATED(matrix)) DEALLOCATE(matrix,dst_add,src_add)
DEALLOCATE( masksrc )
DEALLOCATE(detected_pts)
!
! end loop on vertical levels
!
END DO
ENDIF
CALL Copy_Ncdf_att(TRIM(varname),TRIM(filename),Child_file)
!
......