!---------------------------------------------------------------------- ! NEMO system team, System and Interface for oceanic RElocable Nesting !---------------------------------------------------------------------- ! ! ! PROGRAM: create_coord ! ! DESCRIPTION: !> @brief !> This program create coordinate file. !> !> @details !> Variables are extracted from the input coordinates coarse grid, !> and interpolated to create fine coordinates files. !> !> @author !> J.Paul ! REVISION HISTORY: !> @date Nov, 2013 - Initial Version ! !> @note Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) !> !> @todo !> - add extrapolation (case coordin with mask) !> - add extraction from a grid at fine resolution !---------------------------------------------------------------------- !> @code PROGRAM create_coord ! USE netcdf ! nf90 library USE global ! global variable USE kind ! F90 kind parameter USE logger ! log file manager USE fct ! basic useful function USE date ! date manager USE att ! attribute manager USE dim ! dimension manager USE var ! variable manager USE file ! file manager USE iom ! I/O manager USE dom ! domain manager USE grid ! grid manager USE extrap ! extrapolation manager USE interp ! interpolation manager USE filter ! filter manager USE mpp ! MPP manager USE iom_mpp ! MPP I/O manager IMPLICIT NONE ! local variable CHARACTER(LEN=lc) :: cl_namelist CHARACTER(LEN=lc) :: cl_date INTEGER(i4) :: il_narg INTEGER(i4) :: il_status INTEGER(i4) :: il_fileid INTEGER(i4) :: il_nvar ! INTEGER(i4) , DIMENSION(:,:,:,:) , ALLOCATABLE :: il_value INTEGER(i4) , DIMENSION(ip_maxdim) :: il_rho LOGICAL :: ll_exist TYPE(TATT) :: tl_att TYPE(TDOM) :: tl_dom TYPE(TVAR) , DIMENSION(:) , ALLOCATABLE :: tl_var TYPE(TDIM) , DIMENSION(ip_maxdim) :: tl_dim TYPE(TFILE) :: tl_coord0 TYPE(TFILE) :: tl_fileout TYPE(TMPP) :: tl_mppcoordin ! loop indices INTEGER(i4) :: ji ! namelist variable CHARACTER(LEN=lc) :: cn_logfile = 'create_coord.log' CHARACTER(LEN=lc) :: cn_verbosity = 'warning' CHARACTER(LEN=lc) :: cn_coord0 = '' INTEGER(i4) :: in_perio0 = -1 CHARACTER(LEN=lc) :: cn_varcfg = 'variable.cfg' CHARACTER(LEN=lc), DIMENSION(ig_maxvar) :: cn_varinfo = '' INTEGER(i4) :: in_imin0 = 0 INTEGER(i4) :: in_imax0 = 0 INTEGER(i4) :: in_jmin0 = 0 INTEGER(i4) :: in_jmax0 = 0 INTEGER(i4) :: in_rhoi = 1 INTEGER(i4) :: in_rhoj = 1 CHARACTER(LEN=lc) :: cn_fileout= 'coord_fine.nc' !------------------------------------------------------------------- NAMELIST /namlog/ & !< logger namelist & cn_logfile, & !< log file & cn_verbosity !< logger verbosity NAMELIST /namcfg/ & !< config namelist & cn_varcfg !< variable configuration file NAMELIST /namcrs/ & ! coarse grid namelist & cn_coord0 , & !< coordinate file & in_perio0 !< periodicity index NAMELIST /namvar/ & ! namvar & cn_varinfo !< list of variable and extra information about !< interpolation, extrapolation or filter method to be used. !< (ex: 'votemper/linear/hann/dist_weight','vosaline/cubic' ) NAMELIST /namnst/ & !< nesting namelist & in_imin0, & !< i-direction lower left point indice & in_imax0, & !< i-direction upper right point indice & in_jmin0, & !< j-direction lower left point indice & in_jmax0, & !< j-direction upper right point indice & in_rhoi, & !< refinement factor in i-direction & in_rhoj !< refinement factor in j-direction NAMELIST /namout/ & !< output namelist & cn_fileout !< fine grid coordinate file !------------------------------------------------------------------- !1- namelist !1-1 get namelist il_narg=COMMAND_ARGUMENT_COUNT() !f03 intrinsec IF( il_narg/=1 )THEN PRINT *,"ERROR in create_coord: need a namelist" STOP ELSE CALL GET_COMMAND_ARGUMENT(1,cl_namelist) !f03 intrinsec ENDIF !1-2 read namelist INQUIRE(FILE=TRIM(cl_namelist), EXIST=ll_exist) IF( ll_exist )THEN il_fileid=fct_getunit() OPEN( il_fileid, FILE=TRIM(cl_namelist), & & FORM='FORMATTED', & & ACCESS='SEQUENTIAL', & & STATUS='OLD', & & ACTION='READ', & & IOSTAT=il_status) CALL fct_err(il_status) IF( il_status /= 0 )THEN PRINT *,"ERROR in create_coord: error opening "//TRIM(cl_namelist) STOP ENDIF READ( il_fileid, NML = namlog ) !1-2-1 define log file CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity)) CALL logger_header() READ( il_fileid, NML = namcfg ) !1-2-2 get variable extra information on configuration file CALL var_def_extra(TRIM(cn_varcfg)) READ( il_fileid, NML = namcrs ) READ( il_fileid, NML = namvar ) !1-2-3 add user change in extra information CALL var_chg_extra( cn_varinfo ) READ( il_fileid, NML = namnst ) READ( il_fileid, NML = namout ) CLOSE( il_fileid, IOSTAT=il_status ) CALL fct_err(il_status) IF( il_status /= 0 )THEN CALL logger_error("CREATE COORD: closing "//TRIM(cl_namelist)) ENDIF ELSE PRINT *,"ERROR in create_coord: can't find "//TRIM(cl_namelist) ENDIF !2- open files IF( cn_coord0 /= '' )THEN tl_coord0=file_init(TRIM(cn_coord0),id_perio=in_perio0) CALL iom_open(tl_coord0) ELSE CALL logger_fatal("CREATE COORD: no coarse grid coordinate found. "//& & "check namelist") ENDIF !3- check !3-1 check output file do not already exist INQUIRE(FILE=TRIM(cn_fileout), EXIST=ll_exist) IF( ll_exist )THEN CALL logger_fatal("CREATE COORD: output file "//TRIM(cn_fileout)//& & " already exist.") ENDIF !3-2 check namelist IF( in_imin0 < 1 .OR. in_imax0 < 1 .OR. in_jmin0 < 1 .OR. in_jmax0 < 1)THEN CALL logger_error("CREATE COORD: invalid point indice."//& & " check namelist "//TRIM(cl_namelist)) ENDIF il_rho(:)=1 IF( in_rhoi < 1 .OR. in_rhoj < 1 )THEN CALL logger_error("CREATE COORD: invalid refinement factor."//& & " check namelist "//TRIM(cl_namelist)) ELSE il_rho(jp_I)=in_rhoi il_rho(jp_J)=in_rhoj ENDIF !3-3 check domain validity CALL grid_check_dom(tl_coord0, in_imin0, in_imax0, in_jmin0, in_jmax0 ) !4- compute domain tl_dom=dom_init( tl_coord0, & & in_imin0, in_imax0,& & in_jmin0, in_jmax0 ) ! close file CALL iom_close(tl_coord0) !4-1 add extra band (if possible) to compute interpolation CALL dom_add_extra(tl_dom) !5- read variables on domain (ugly way to do it, have to work on it) !5-1 init mpp structure tl_mppcoordin=mpp_init(tl_coord0) CALL file_clean(tl_coord0) !5-2 get processor to be used CALL mpp_get_use( tl_mppcoordin, tl_dom ) !5-3 open mpp files CALL iom_mpp_open(tl_mppcoordin) !5-4 fill variable value on domain CALL iom_mpp_fill_var(tl_mppcoordin, tl_dom) !5-5 close mpp files CALL iom_mpp_close(tl_mppcoordin) il_nvar=tl_mppcoordin%t_proc(1)%i_nvar ALLOCATE( tl_var(il_nvar) ) DO ji=1,il_nvar tl_var(ji)=tl_mppcoordin%t_proc(1)%t_var(ji) !7- interpolate variables CALL create_coord_interp( tl_var(ji), il_rho(:) ) !6- remove extraband added to domain CALL dom_del_extra( tl_var(ji), tl_dom, il_rho(:) ) !7- add ghost cell CALL grid_add_ghost(tl_var(ji),tl_dom%i_ighost,tl_dom%i_jghost) !8- filter CALL filter_fill_value(tl_var(ji)) ENDDO !9- clean DO ji=1,il_nvar CALL var_clean(tl_mppcoordin%t_proc(1)%t_var(ji)) ENDDO CALL mpp_clean(tl_mppcoordin) !10- create file tl_fileout=file_init(TRIM(cn_fileout)) !10-1 add dimension ! save biggest dimension tl_dim(:)=var_max_dim(tl_var(:)) DO ji=1,ip_maxdim IF( tl_dim(ji)%l_use ) CALL file_add_dim(tl_fileout, tl_dim(ji)) ENDDO !10-2 add variables DO ji=1,il_nvar CALL file_add_var(tl_fileout, tl_var(ji)) ENDDO !10-3 add some attribute tl_att=att_init("Created_by","SIREN create_coord") CALL file_add_att(tl_fileout, tl_att) cl_date=date_print(date_now()) tl_att=att_init("Creation_date",cl_date) CALL file_add_att(tl_fileout, tl_att) tl_att=att_init("source_file",TRIM(fct_basename(cn_coord0))) CALL file_add_att(tl_fileout, tl_att) tl_att=att_init("source_i-indices",(/in_imin0,in_imax0/)) CALL file_add_att(tl_fileout, tl_att) tl_att=att_init("source_j-indices",(/in_jmin0,in_jmax0/)) CALL file_add_att(tl_fileout, tl_att) !10-4 create file CALL iom_create(tl_fileout) !10-5 write file CALL iom_write_file(tl_fileout) !10-6 close file CALL iom_close(tl_fileout) !11- clean DO ji=1,il_nvar CALL var_clean(tl_var(ji)) ENDDO CALL file_clean(tl_fileout) DEALLOCATE( tl_var) ! close log file CALL logger_footer() CALL logger_close() !> @endcode CONTAINS !------------------------------------------------------------------- !> @brief !> This subroutine !> !> @details !> !> @author J.Paul !> - Nov, 2013- Initial Version !> !> @param[in] !> @todo !------------------------------------------------------------------- !> @code SUBROUTINE create_coord_interp( td_var, & & id_rho, & & id_iext, id_jext) IMPLICIT NONE ! Argument TYPE(TVAR) , INTENT(INOUT) :: td_var INTEGER(i4), DIMENSION(:), INTENT(IN ) :: id_rho INTEGER(i4), INTENT(IN ), OPTIONAL :: id_iext INTEGER(i4), INTENT(IN ), OPTIONAL :: id_jext ! local variable TYPE(TVAR) :: tl_mask TYPE(TVAR) :: tl_var INTEGER(i1), DIMENSION(:,:,:,:), ALLOCATABLE :: bl_mask INTEGER(i4), DIMENSION(2,2) :: il_offset INTEGER(i4) :: il_iext INTEGER(i4) :: il_jext ! loop indices !---------------------------------------------------------------- ! copy variable tl_var=td_var !WARNING: two extrabands are required for cubic interpolation il_iext=2 IF( PRESENT(id_iext) ) il_iext=id_iext il_jext=2 IF( PRESENT(id_jext) ) il_jext=id_jext IF( il_iext < 2 .AND. td_var%c_interp(1) == 'cubic' )THEN CALL logger_warn("CREATE COORD INTERP: at least extrapolation "//& & "on two points are required with cubic interpolation ") il_iext=2 ENDIF IF( il_jext < 2 .AND. td_var%c_interp(1) == 'cubic' )THEN CALL logger_warn("CREATE COORD INTERP: at least extrapolation "//& & "on two points are required with cubic interpolation ") il_jext=2 ENDIF !1- work on mask !1-1 create mask ALLOCATE(bl_mask(tl_var%t_dim(1)%i_len, & & tl_var%t_dim(2)%i_len, & & tl_var%t_dim(3)%i_len, & & tl_var%t_dim(4)%i_len) ) bl_mask(:,:,:,:)=1 WHERE(tl_var%d_value(:,:,:,:)==tl_var%d_fill) bl_mask(:,:,:,:)=0 SELECT CASE(TRIM(tl_var%c_point)) CASE DEFAULT ! 'T' tl_mask=var_init('tmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) CASE('U') tl_mask=var_init('umask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) CASE('V') tl_mask=var_init('vmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) CASE('F') tl_mask=var_init('fmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) END SELECT DEALLOCATE(bl_mask) !1-2 interpolate mask il_offset(:,:)=1 CALL interp_fill_value( tl_mask, id_rho(:), & & id_offset=il_offset(:,:) ) !2- work on variable !2-0 add extraband CALL extrap_add_extrabands(tl_var, il_iext, il_jext) !2-1 extrapolate variable CALL extrap_fill_value( tl_var, id_iext=il_iext, id_jext=il_jext ) !2-2 interpolate variable il_offset(:,:)=1 CALL interp_fill_value( tl_var, id_rho(:), & & id_offset=il_offset(:,:)) !2-3 remove extraband CALL extrap_del_extrabands(tl_var, il_iext*id_rho(jp_I), il_jext*id_rho(jp_J)) !3- keep original mask WHERE( tl_mask%d_value(:,:,:,:) == 0 ) tl_var%d_value(:,:,:,:)=tl_var%d_fill END WHERE !4- save result td_var=tl_var ! clean variable structure CALL var_clean(tl_mask) CALL var_clean(tl_var) END SUBROUTINE create_coord_interp !> @endcode END PROGRAM create_coord