Ignore:
Timestamp:
2015-07-15T17:46:12+02:00 (5 years ago)
Author:
andrewryan
Message:

merged in latest version of trunk alongside changes to SAO_SRC to be compatible with latest OBS

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/TOOLS/SIREN/src/create_coord.f90

    r4213 r5600  
    77! 
    88! DESCRIPTION: 
     9!> @file 
    910!> @brief  
    10 !> This program create coordinate file. 
     11!> This program create fine grid coordinate file. 
    1112!> 
    1213!> @details 
    13 !> Variables are extracted from the input coordinates coarse grid,  
    14 !> and interpolated to create fine coordinates files. 
    15 !> 
    16 !> @author 
    17 !> J.Paul 
     14!> @section sec1 method 
     15!>    All variables from the input coordinates coarse grid file, are extracted 
     16!>    and interpolated to create fine grid coordinates files.<br/> 
     17!>    @note  
     18!>       interpolation method could be different for each variable. 
     19!> 
     20!> @section sec2 how to 
     21!>    to create fine grid coordinates files:<br/> 
     22!> @code{.sh} 
     23!>    ./SIREN/bin/create_coord create_coord.nam 
     24!> @endcode 
     25!>     
     26!>    create_coord.nam comprise 6 namelists:<br/> 
     27!>       - logger namelist (namlog) 
     28!>       - config namelist (namcfg) 
     29!>       - coarse grid namelist (namcrs) 
     30!>       - variable namelist (namvar) 
     31!>       - nesting namelist (namnst) 
     32!>       - output namelist (namout) 
     33!>     
     34!>    @note  
     35!>       All namelists have to be in file create_coord.nam,  
     36!>       however variables of those namelists are all optional. 
     37!> 
     38!>    * _logger namelist (namlog)_:<br/> 
     39!>       - cn_logfile   : log filename 
     40!>       - cn_verbosity : verbosity ('trace','debug','info', 
     41!> 'warning','error','fatal') 
     42!>       - in_maxerror  : maximum number of error allowed 
     43!> 
     44!>    * _config namelist (namcfg)_:<br/> 
     45!>       - cn_varcfg : variable configuration file  
     46!> (see ./SIREN/cfg/variable.cfg) 
     47!> 
     48!>    * _coarse grid namelist (namcrs)_:<br/> 
     49!>       - cn_coord0 : coordinate file 
     50!>       - in_perio0 : NEMO periodicity index (see Model Boundary Condition in 
     51!> [NEMO documentation](http://www.nemo-ocean.eu/About-NEMO/Reference-manuals)) 
     52!> 
     53!>    * _variable namelist (namvar)_:<br/> 
     54!>       - cn_varinfo : list of variable and extra information about request(s) 
     55!> to be used.<br/> 
     56!>          each elements of *cn_varinfo* is a string character.<br/> 
     57!>          it is composed of the variable name follow by ':',  
     58!>          then request(s) to be used on this variable.<br/>  
     59!>          request could be: 
     60!>             - interpolation method 
     61!>             - extrapolation method 
     62!>             - filter method 
     63!>  
     64!>                requests must be separated by ';' .<br/> 
     65!>                order of requests does not matter.<br/> 
     66!> 
     67!>          informations about available method could be find in @ref interp, 
     68!>          @ref extrap and @ref filter modules.<br/> 
     69!> 
     70!>          Example: 'votemper: linear; hann(2,3); dist_weight',  
     71!>          'vosaline: cubic'<br/> 
     72!>          @note  
     73!>             If you do not specify a method which is required,  
     74!>             default one is applied. 
     75!> 
     76!>    * _nesting namelist (namnst)_:<br/> 
     77!>       - in_imin0 : i-direction lower left  point indice  
     78!> of coarse grid subdomain to be used 
     79!>       - in_imax0 : i-direction upper right point indice 
     80!> of coarse grid subdomain to be used 
     81!>       - in_jmin0 : j-direction lower left  point indice 
     82!> of coarse grid subdomain to be used 
     83!>       - in_jmax0 : j-direction upper right point indice 
     84!> of coarse grid subdomain to be used 
     85!>       - in_rhoi  : refinement factor in i-direction 
     86!>       - in_rhoj  : refinement factor in j-direction<br/> 
     87!> 
     88!>       \image html  grid_zoom_40.png  
     89!>       \image latex grid_zoom_40.png 
     90!> 
     91!>    * _output namelist (namout)_: 
     92!>       - cn_fileout : output coordinate file 
     93!> 
     94!> @author J.Paul 
    1895! REVISION HISTORY: 
    19 !> @date Nov, 2013 - Initial Version 
    20 ! 
     96!> @date November, 2013 - Initial Version 
     97!> @date September, 2014 
     98!> - add header for user 
     99!> - compute offset considering grid point 
     100!> - add global attributes in output file 
     101!> 
    21102!> @note Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    22 !> 
    23 !> @todo 
    24 !> - add extrapolation (case coordin with mask) 
    25 !> - add extraction from a grid at fine resolution 
    26103!---------------------------------------------------------------------- 
    27 !> @code 
    28104PROGRAM create_coord 
    29105 
    30 !   USE netcdf                          ! nf90 library 
    31106   USE global                          ! global variable 
    32107   USE kind                            ! F90 kind parameter 
     
    39114   USE file                            ! file manager 
    40115   USE iom                             ! I/O manager 
    41    USE dom                             ! domain manager 
    42116   USE grid                            ! grid manager 
    43117   USE extrap                          ! extrapolation manager 
     
    45119   USE filter                          ! filter manager 
    46120   USE mpp                             ! MPP manager 
     121   USE dom                             ! domain manager 
    47122   USE iom_mpp                         ! MPP I/O manager 
     123   USE iom_dom                         ! DOM I/O manager 
    48124 
    49125   IMPLICIT NONE 
     
    56132   INTEGER(i4)                                          :: il_status 
    57133   INTEGER(i4)                                          :: il_fileid 
     134   INTEGER(i4)                                          :: il_attid 
     135   INTEGER(i4)                                          :: il_ind 
    58136   INTEGER(i4)                                          :: il_nvar 
    59 !   INTEGER(i4)      , DIMENSION(:,:,:,:)  , ALLOCATABLE :: il_value 
     137   INTEGER(i4)                                          :: il_ew 
    60138   INTEGER(i4)      , DIMENSION(ip_maxdim)              :: il_rho 
    61  
     139   INTEGER(i4)      , DIMENSION(2,2,ip_npoint)          :: il_offset 
    62140 
    63141   LOGICAL                                              :: ll_exist 
     
    71149   TYPE(TDIM)       , DIMENSION(ip_maxdim)              :: tl_dim 
    72150 
    73    TYPE(TFILE)                                          :: tl_coord0 
     151   TYPE(TMPP)                                           :: tl_coord0 
    74152   TYPE(TFILE)                                          :: tl_fileout 
    75153 
    76    TYPE(TMPP)                                           :: tl_mppcoordin 
     154   ! check  
     155!   INTEGER(i4)                                          :: il_imin0 
     156!   INTEGER(i4)                                          :: il_imax0 
     157!   INTEGER(i4)                                          :: il_jmin0 
     158!   INTEGER(i4)                                          :: il_jmax0 
     159!   INTEGER(i4)      , DIMENSION(2,2)                    :: il_ind2 
     160!   TYPE(TMPP)                                           :: tl_mppout 
    77161 
    78162   ! loop indices 
    79163   INTEGER(i4) :: ji 
     164   INTEGER(i4) :: jj 
    80165 
    81166   ! namelist variable 
    82167   CHARACTER(LEN=lc) :: cn_logfile = 'create_coord.log'  
    83168   CHARACTER(LEN=lc) :: cn_verbosity = 'warning'  
     169   INTEGER(i4)       :: in_maxerror = 5 
    84170 
    85171   CHARACTER(LEN=lc) :: cn_coord0 = ''  
    86172   INTEGER(i4)       :: in_perio0 = -1 
    87173 
    88    CHARACTER(LEN=lc) :: cn_varcfg = 'variable.cfg'  
    89  
    90    CHARACTER(LEN=lc), DIMENSION(ig_maxvar) :: cn_varinfo = '' 
     174   CHARACTER(LEN=lc) :: cn_varcfg = '../cfg/variable.cfg'  
     175 
     176   CHARACTER(LEN=lc), DIMENSION(ip_maxvar) :: cn_varinfo = '' 
    91177 
    92178   INTEGER(i4)       :: in_imin0 = 0 
     
    100186   !------------------------------------------------------------------- 
    101187 
    102    NAMELIST /namlog/ & !< logger namelist 
    103    &  cn_logfile,    & !< log file 
    104    &  cn_verbosity     !< logger verbosity 
    105  
    106    NAMELIST /namcfg/ &  !< config namelist 
     188   NAMELIST /namlog/ &  !  logger namelist 
     189   &  cn_logfile,    &  !< logger file name 
     190   &  cn_verbosity,  &  !< logger verbosity 
     191   &  in_maxerror       !< logger maximum error 
     192 
     193   NAMELIST /namcfg/ &  !  config namelist 
    107194   &  cn_varcfg         !< variable configuration file 
    108195 
    109    NAMELIST /namcrs/ &  ! coarse grid namelist 
     196   NAMELIST /namcrs/ &  !  coarse grid namelist 
    110197   &  cn_coord0 , &     !< coordinate file 
    111198   &  in_perio0         !< periodicity index 
    112199 
    113    NAMELIST /namvar/ &  ! namvar 
     200   NAMELIST /namvar/ &  !  variable namelist 
    114201   &  cn_varinfo        !< list of variable and extra information about  
    115202                        !< interpolation, extrapolation or filter method to be used.  
    116                         !< (ex: 'votemper/linear/hann/dist_weight','vosaline/cubic' )  
     203                        !< (ex: 'votemper:linear,hann,dist_weight','vosaline:cubic' )  
    117204    
    118    NAMELIST /namnst/ &  !< nesting namelist 
     205   NAMELIST /namnst/ &  !  nesting namelist 
    119206   &  in_imin0,   &     !< i-direction lower left  point indice  
    120207   &  in_imax0,   &     !< i-direction upper right point indice 
     
    124211   &  in_rhoj           !< refinement factor in j-direction 
    125212 
    126    NAMELIST /namout/ &  !< output namelist 
    127    &  cn_fileout       !< fine grid coordinate file    
     213   NAMELIST /namout/ &  !  output namelist 
     214   &  cn_fileout        !< fine grid coordinate file    
    128215   !------------------------------------------------------------------- 
    129216 
    130    !1- namelist 
    131    !1-1 get namelist 
     217   ! namelist 
     218   ! get namelist 
    132219   il_narg=COMMAND_ARGUMENT_COUNT() !f03 intrinsec 
    133220   IF( il_narg/=1 )THEN 
     
    138225   ENDIF 
    139226    
    140    !1-2 read namelist 
     227   ! read namelist 
    141228   INQUIRE(FILE=TRIM(cl_namelist), EXIST=ll_exist) 
    142229   IF( ll_exist )THEN 
    143        
     230  
    144231      il_fileid=fct_getunit() 
    145232 
     
    157244 
    158245      READ( il_fileid, NML = namlog ) 
    159       !1-2-1 define log file 
    160       CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity)) 
     246      ! define logger file 
     247      CALL logger_open(TRIM(cn_logfile),TRIM(cn_verbosity),in_maxerror) 
    161248      CALL logger_header() 
    162249 
    163250      READ( il_fileid, NML = namcfg ) 
    164       !1-2-2 get variable extra information on configuration file 
     251      ! get variable extra information on configuration file 
    165252      CALL var_def_extra(TRIM(cn_varcfg)) 
    166253 
    167254      READ( il_fileid, NML = namcrs ) 
    168255      READ( il_fileid, NML = namvar ) 
    169       !1-2-3 add user change in extra information 
     256      ! add user change in extra information 
    170257      CALL var_chg_extra( cn_varinfo ) 
    171258 
     
    182269 
    183270      PRINT *,"ERROR in create_coord: can't find "//TRIM(cl_namelist) 
    184  
    185    ENDIF 
    186  
    187    !2- open files 
     271      STOP 
     272 
     273   ENDIF 
     274 
     275   ! open files 
    188276   IF( cn_coord0 /= '' )THEN 
    189       tl_coord0=file_init(TRIM(cn_coord0),id_perio=in_perio0) 
    190       CALL iom_open(tl_coord0) 
     277      tl_coord0=mpp_init( file_init(TRIM(cn_coord0)), id_perio=in_perio0) 
     278      CALL grid_get_info(tl_coord0) 
    191279   ELSE 
    192        CALL logger_fatal("CREATE COORD: no coarse grid coordinate found. "//& 
     280      CALL logger_fatal("CREATE COORD: no coarse grid coordinate found. "//& 
    193281      &     "check namelist")       
    194282   ENDIF 
    195283 
    196    !3- check 
    197    !3-1 check output file do not already exist 
     284   ! check 
     285   ! check output file do not already exist 
    198286   INQUIRE(FILE=TRIM(cn_fileout), EXIST=ll_exist) 
    199287   IF( ll_exist )THEN 
     
    202290   ENDIF 
    203291 
    204    !3-2 check namelist 
    205    IF( in_imin0 < 1 .OR. in_imax0 < 1 .OR. in_jmin0 < 1 .OR. in_jmax0 < 1)THEN 
    206       CALL logger_error("CREATE COORD: invalid point indice."//& 
     292   ! check nesting parameters 
     293   IF( in_imin0 < 0 .OR. in_imax0 < 0 .OR. in_jmin0 < 0 .OR. in_jmax0 < 0)THEN 
     294      CALL logger_fatal("CREATE COORD: invalid points indices."//& 
    207295      &  " check namelist "//TRIM(cl_namelist)) 
    208296   ENDIF 
     
    215303      il_rho(jp_I)=in_rhoi 
    216304      il_rho(jp_J)=in_rhoj       
    217    ENDIF 
    218  
    219    !3-3 check domain validity 
     305 
     306      il_offset(:,:,:)=create_coord_get_offset(il_rho(:)) 
     307 
     308   ENDIF 
     309 
     310   ! check domain validity 
    220311   CALL grid_check_dom(tl_coord0, in_imin0, in_imax0, in_jmin0, in_jmax0 ) 
    221312 
    222    !4- compute domain 
     313   ! compute domain 
    223314   tl_dom=dom_init( tl_coord0,         & 
    224315   &                in_imin0, in_imax0,& 
    225316   &                in_jmin0, in_jmax0 ) 
    226317 
    227    ! close file 
    228    CALL iom_close(tl_coord0) 
    229  
    230    !4-1 add extra band (if possible) to compute interpolation 
     318   ! add extra band (if need be) to compute interpolation 
    231319   CALL dom_add_extra(tl_dom) 
    232320 
    233    !5- read variables on domain (ugly way to do it, have to work on it) 
    234    !5-1 init mpp structure 
    235    tl_mppcoordin=mpp_init(tl_coord0) 
    236  
    237    CALL file_clean(tl_coord0) 
    238  
    239    !5-2 get processor to be used 
    240    CALL mpp_get_use( tl_mppcoordin, tl_dom ) 
    241  
    242    !5-3 open mpp files 
    243    CALL iom_mpp_open(tl_mppcoordin) 
    244  
    245    !5-4 fill variable value on domain 
    246    CALL iom_mpp_fill_var(tl_mppcoordin, tl_dom) 
    247  
    248    !5-5 close mpp files 
    249    CALL iom_mpp_close(tl_mppcoordin) 
    250  
    251    il_nvar=tl_mppcoordin%t_proc(1)%i_nvar 
     321   ! open mpp files 
     322   CALL iom_dom_open(tl_coord0, tl_dom) 
     323 
     324   il_nvar=tl_coord0%t_proc(1)%i_nvar 
    252325   ALLOCATE( tl_var(il_nvar) ) 
    253326   DO ji=1,il_nvar 
    254327 
    255       tl_var(ji)=tl_mppcoordin%t_proc(1)%t_var(ji) 
    256       !7- interpolate variables 
    257       CALL create_coord_interp( tl_var(ji), il_rho(:) ) 
    258  
    259       !6- remove extraband added to domain 
    260       CALL dom_del_extra( tl_var(ji), tl_dom, il_rho(:) ) 
    261  
    262       !7- add ghost cell 
    263       CALL grid_add_ghost(tl_var(ji),tl_dom%i_ighost,tl_dom%i_jghost)       
    264  
    265       !8- filter 
     328      tl_var(ji)=iom_dom_read_var(tl_coord0, & 
     329      &                          TRIM(tl_coord0%t_proc(1)%t_var(ji)%c_name),& 
     330      &                          tl_dom) 
     331 
     332      SELECT CASE(TRIM(tl_var(ji)%c_point)) 
     333         CASE('T') 
     334            jj=jp_T 
     335         CASE('U') 
     336            jj=jp_U 
     337         CASE('V') 
     338            jj=jp_V 
     339         CASE('F') 
     340            jj=jp_F 
     341      END SELECT 
     342 
     343      ! interpolate variables 
     344      CALL create_coord_interp( tl_var(ji), il_rho(:), & 
     345      &                         il_offset(:,:,jj) ) 
     346 
     347      ! remove extraband added to domain 
     348      CALL dom_del_extra( tl_var(ji), tl_dom, il_rho(:), .true. ) 
     349 
     350      ! do not add ghost cell.  
     351      ! ghost cell already replace by value for coordinates  
     352      ! CALL grid_add_ghost(tl_var(ji),tl_dom%i_ghost(:,:)) 
     353 
     354      ! filter 
    266355      CALL filter_fill_value(tl_var(ji))       
    267356 
    268357   ENDDO 
    269358 
    270    !9- clean 
    271    DO ji=1,il_nvar 
    272       CALL var_clean(tl_mppcoordin%t_proc(1)%t_var(ji)) 
    273    ENDDO 
    274    CALL mpp_clean(tl_mppcoordin) 
    275  
    276    !10- create file 
     359   ! close mpp files 
     360   CALL iom_dom_close(tl_coord0) 
     361 
     362   ! clean 
     363   CALL mpp_clean(tl_coord0) 
     364 
     365   ! create file 
    277366   tl_fileout=file_init(TRIM(cn_fileout)) 
    278367 
    279    !10-1 add dimension 
     368   ! add dimension 
    280369   ! save biggest dimension 
    281370   tl_dim(:)=var_max_dim(tl_var(:)) 
     
    285374   ENDDO 
    286375 
    287    !10-2 add variables 
    288  
     376   ! add variables 
    289377   DO ji=1,il_nvar 
    290378      CALL file_add_var(tl_fileout, tl_var(ji)) 
    291379   ENDDO 
    292380 
    293    !10-3 add some attribute 
     381   ! recompute some attribute 
     382 
     383   ! add some attribute 
    294384   tl_att=att_init("Created_by","SIREN create_coord") 
    295385   CALL file_add_att(tl_fileout, tl_att) 
     
    299389   CALL file_add_att(tl_fileout, tl_att) 
    300390 
    301    tl_att=att_init("source_file",TRIM(fct_basename(cn_coord0))) 
     391   tl_att=att_init("src_file",TRIM(fct_basename(cn_coord0))) 
    302392   CALL file_add_att(tl_fileout, tl_att)    
    303393 
    304    tl_att=att_init("source_i-indices",(/in_imin0,in_imax0/)) 
     394   tl_att=att_init("src_i_indices",(/in_imin0,in_imax0/)) 
    305395   CALL file_add_att(tl_fileout, tl_att)    
    306    tl_att=att_init("source_j-indices",(/in_jmin0,in_jmax0/)) 
    307    CALL file_add_att(tl_fileout, tl_att)    
    308  
    309    !10-4 create file 
     396   tl_att=att_init("src_j_indices",(/in_jmin0,in_jmax0/)) 
     397   CALL file_add_att(tl_fileout, tl_att) 
     398   IF( .NOT. ALL(il_rho(:)==1) )THEN 
     399      tl_att=att_init("refinment_factor",(/il_rho(jp_I),il_rho(jp_J)/)) 
     400      CALL file_add_att(tl_fileout, tl_att) 
     401   ENDIF 
     402 
     403   ! add attribute periodicity 
     404   il_attid=0 
     405   IF( ASSOCIATED(tl_fileout%t_att) )THEN 
     406      il_attid=att_get_index(tl_fileout%t_att(:),'periodicity') 
     407   ENDIF 
     408   IF( tl_dom%i_perio >= 0 .AND. il_attid == 0 )THEN 
     409      tl_att=att_init('periodicity',tl_dom%i_perio) 
     410      CALL file_add_att(tl_fileout,tl_att) 
     411   ENDIF 
     412 
     413   ! add attribute east west overlap 
     414   il_attid=0 
     415   IF( ASSOCIATED(tl_fileout%t_att) )THEN 
     416      il_attid=att_get_index(tl_fileout%t_att(:),'ew_overlap') 
     417   ENDIF 
     418   IF( il_attid == 0 )THEN 
     419      il_ind=var_get_index(tl_fileout%t_var(:),'longitude') 
     420      il_ew=grid_get_ew_overlap(tl_fileout%t_var(il_ind)) 
     421      IF( il_ew >= 0 )THEN 
     422         tl_att=att_init('ew_overlap',il_ew) 
     423         CALL file_add_att(tl_fileout,tl_att) 
     424      ENDIF 
     425   ENDIF 
     426 
     427   ! create file 
    310428   CALL iom_create(tl_fileout) 
    311429 
    312    !10-5 write file 
     430   ! write file 
    313431   CALL iom_write_file(tl_fileout) 
    314432 
    315    !10-6 close file 
     433   ! close file 
    316434   CALL iom_close(tl_fileout) 
    317435 
    318    !11- clean 
    319    DO ji=1,il_nvar 
    320       CALL var_clean(tl_var(ji)) 
    321    ENDDO 
     436   ! clean 
     437   CALL att_clean(tl_att) 
     438   CALL var_clean(tl_var(:)) 
     439   DEALLOCATE( tl_var)  
     440 
    322441   CALL file_clean(tl_fileout) 
    323442 
    324    DEALLOCATE( tl_var)  
     443!   ! check domain 
     444!   tl_coord0=mpp_init( file_init(TRIM(cn_coord0)), id_perio=in_perio0) 
     445!   tl_mppout=mpp_init( file_init(TRIM(cn_fileout)) ) 
     446!   CALL grid_get_info(tl_coord0) 
     447!   CALL iom_mpp_open(tl_mppout) 
     448! 
     449!   il_ind2(:,:)=grid_get_coarse_index( tl_coord0, tl_mppout, & 
     450!   &                                   id_rho=il_rho(:) ) 
     451! 
     452!   il_imin0=il_ind2(1,1) ; il_imax0=il_ind2(1,2) 
     453!   il_jmin0=il_ind2(2,1) ; il_jmax0=il_ind2(2,2) 
     454! 
     455!   IF( il_imin0 /= in_imin0 .OR. & 
     456!   &   il_imax0 /= in_imax0 .OR. & 
     457!   &   il_jmin0 /= in_jmin0 .OR. & 
     458!   &   il_jmax0 /= in_jmax0 )THEN 
     459!      CALL logger_debug("CREATE COORD: output indices ("//& 
     460!      &                 TRIM(fct_str(il_imin0))//","//& 
     461!      &                 TRIM(fct_str(il_imax0))//") ("//& 
     462!      &                 TRIM(fct_str(il_jmin0))//","//& 
     463!      &                 TRIM(fct_str(il_jmax0))//")" )  
     464!      CALL logger_debug("CREATE COORD: input indices ("//& 
     465!      &                 TRIM(fct_str(in_imin0))//","//& 
     466!      &                 TRIM(fct_str(in_imax0))//") ("//& 
     467!      &                 TRIM(fct_str(in_jmin0))//","//& 
     468!      &                 TRIM(fct_str(in_jmax0))//")" )  
     469!      CALL logger_fatal("CREATE COORD: output domain not confrom "//& 
     470!      &                 "with input indices") 
     471!   ENDIF 
     472! 
     473!   CALL iom_mpp_close(tl_coord0) 
     474!   CALL iom_mpp_close(tl_mppout) 
    325475 
    326476   ! close log file 
    327477   CALL logger_footer() 
    328    CALL logger_close()    
    329  
    330 !> @endcode 
     478   CALL logger_close()  
     479 
    331480CONTAINS 
    332481   !------------------------------------------------------------------- 
    333482   !> @brief 
    334    !> This subroutine 
     483   !> This function compute offset over Arakawa grid points,  
     484   !> given refinement factor. 
     485   !>  
     486   !> @author J.Paul 
     487   !> @date August, 2014 - Initial Version 
     488   !> 
     489   !> @param[in] id_rho array of refinement factor 
     490   !> @return array of offset 
     491   !------------------------------------------------------------------- 
     492   FUNCTION create_coord_get_offset( id_rho ) 
     493      IMPLICIT NONE 
     494      ! Argument       
     495      INTEGER(i4), DIMENSION(:), INTENT(IN) :: id_rho 
     496 
     497      ! function 
     498      INTEGER(i4), DIMENSION(2,2,ip_npoint) :: create_coord_get_offset 
     499      ! local variable 
     500      ! loop indices 
     501      !---------------------------------------------------------------- 
     502 
     503      ! case 'T' 
     504      create_coord_get_offset(jp_I,:,jp_T)=FLOOR(REAL(id_rho(jp_I)-1,dp)*0.5) 
     505      create_coord_get_offset(jp_J,:,jp_T)=FLOOR(REAL(id_rho(jp_J)-1,dp)*0.5) 
     506      ! case 'U' 
     507      create_coord_get_offset(jp_I,1,jp_U)=0 
     508      create_coord_get_offset(jp_I,2,jp_U)=id_rho(jp_I)-1 
     509      create_coord_get_offset(jp_J,:,jp_U)=FLOOR(REAL(id_rho(jp_J)-1,dp)*0.5) 
     510      ! case 'V' 
     511      create_coord_get_offset(jp_I,:,jp_V)=FLOOR(REAL(id_rho(jp_I)-1,dp)*0.5) 
     512      create_coord_get_offset(jp_J,1,jp_V)=0 
     513      create_coord_get_offset(jp_J,2,jp_V)=id_rho(jp_J)-1 
     514      ! case 'F' 
     515      create_coord_get_offset(jp_I,1,jp_F)=0 
     516      create_coord_get_offset(jp_I,2,jp_F)=id_rho(jp_I)-1 
     517      create_coord_get_offset(jp_J,1,jp_F)=0 
     518      create_coord_get_offset(jp_J,2,jp_F)=id_rho(jp_J)-1 
     519 
     520 
     521   END FUNCTION create_coord_get_offset 
     522   !------------------------------------------------------------------- 
     523   !> @brief 
     524   !> This subroutine interpolate variable, given refinment factor. 
    335525   !>  
    336526   !> @details  
     527   !>  Optionaly, you could specify number of points  
     528   !>    to be extrapolated in i- and j-direction.<br/> 
     529   !>  variable mask is first computed (using _FillValue) and interpolated.<br/> 
     530   !>  variable is then extrapolated, and interpolated.<br/>  
     531   !>  Finally interpolated mask is applied on refined variable. 
    337532   !> 
    338533   !> @author J.Paul 
    339    !> - Nov, 2013- Initial Version 
     534   !> @date November, 2013 - Initial Version 
    340535   !> 
    341    !> @param[in]  
    342    !> @todo  
     536   !> @param[inout] td_var variable strcuture  
     537   !> @param[in] id_rho    array of refinement factor 
     538   !> @param[in] id_offset offset between fine grid and coarse grid 
     539   !> @param[in] id_iext   number of points to be extrapolated in i-direction 
     540   !> @param[in] id_jext   number of points to be extrapolated in j-direction 
    343541   !------------------------------------------------------------------- 
    344    !> @code 
    345542   SUBROUTINE create_coord_interp( td_var,          & 
    346543   &                               id_rho,          & 
     544   &                               id_offset,       & 
    347545   &                               id_iext, id_jext) 
    348546 
     
    350548 
    351549      ! Argument 
    352       TYPE(TVAR) ,               INTENT(INOUT) :: td_var 
    353       INTEGER(i4), DIMENSION(:), INTENT(IN   ) :: id_rho 
    354       INTEGER(i4),               INTENT(IN   ), OPTIONAL :: id_iext 
    355       INTEGER(i4),               INTENT(IN   ), OPTIONAL :: id_jext 
     550      TYPE(TVAR) ,                 INTENT(INOUT) :: td_var 
     551      INTEGER(i4), DIMENSION(:)  , INTENT(IN   ) :: id_rho 
     552      INTEGER(i4), DIMENSION(:,:), INTENT(IN)    :: id_offset 
     553      INTEGER(i4),                 INTENT(IN   ), OPTIONAL :: id_iext 
     554      INTEGER(i4),                 INTENT(IN   ), OPTIONAL :: id_jext 
    356555 
    357556      ! local variable 
    358557      TYPE(TVAR)  :: tl_mask 
    359       TYPE(TVAR)  :: tl_var 
    360558 
    361559      INTEGER(i1), DIMENSION(:,:,:,:), ALLOCATABLE :: bl_mask 
    362  
    363       INTEGER(i4), DIMENSION(2,2) :: il_offset 
    364560 
    365561      INTEGER(i4) :: il_iext 
     
    369565      !---------------------------------------------------------------- 
    370566 
    371       ! copy variable 
    372       tl_var=td_var 
     567      IF( ANY(SHAPE(id_offset(:,:)) /= 2) )THEN 
     568         CALL logger_error("CREATE COORD INTERP: invalid dimension of "//& 
     569         &                 "offset array") 
     570      ENDIF 
    373571 
    374572      !WARNING: two extrabands are required for cubic interpolation 
     
    391589      ENDIF 
    392590 
    393       !1- work on mask 
    394       !1-1 create mask 
    395       ALLOCATE(bl_mask(tl_var%t_dim(1)%i_len, & 
    396       &                tl_var%t_dim(2)%i_len, & 
    397       &                tl_var%t_dim(3)%i_len, & 
    398       &                tl_var%t_dim(4)%i_len) ) 
    399  
    400       bl_mask(:,:,:,:)=1 
    401       WHERE(tl_var%d_value(:,:,:,:)==tl_var%d_fill) bl_mask(:,:,:,:)=0       
    402  
    403       SELECT CASE(TRIM(tl_var%c_point)) 
    404       CASE DEFAULT ! 'T' 
    405          tl_mask=var_init('tmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) 
    406       CASE('U') 
    407          tl_mask=var_init('umask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) 
    408       CASE('V') 
    409          tl_mask=var_init('vmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) 
    410       CASE('F') 
    411          tl_mask=var_init('fmask', bl_mask(:,:,:,:), td_dim=tl_var%t_dim(:)) 
    412       END SELECT          
    413  
    414       DEALLOCATE(bl_mask) 
    415  
    416       !1-2 interpolate mask 
    417       il_offset(:,:)=1 
    418       CALL interp_fill_value( tl_mask, id_rho(:), & 
    419       &                       id_offset=il_offset(:,:) ) 
    420  
    421       !2- work on variable 
    422       !2-0 add extraband 
    423       CALL extrap_add_extrabands(tl_var, il_iext, il_jext) 
    424  
    425       !2-1 extrapolate variable 
    426       CALL extrap_fill_value( tl_var, id_iext=il_iext, id_jext=il_jext ) 
    427  
    428       !2-2 interpolate variable 
    429       il_offset(:,:)=1 
    430       CALL interp_fill_value( tl_var, id_rho(:), & 
    431       &                       id_offset=il_offset(:,:)) 
    432  
    433       !2-3 remove extraband 
    434       CALL extrap_del_extrabands(tl_var, il_iext*id_rho(jp_I), il_jext*id_rho(jp_J)) 
    435  
    436       !3- keep original mask  
    437       WHERE( tl_mask%d_value(:,:,:,:) == 0 ) 
    438          tl_var%d_value(:,:,:,:)=tl_var%d_fill 
    439       END WHERE 
    440  
    441       !4- save result 
    442       td_var=tl_var 
     591      IF( ANY(id_rho(:)>1) )THEN 
     592         ! work on mask 
     593         ! create mask 
     594         ALLOCATE(bl_mask(td_var%t_dim(1)%i_len, & 
     595         &                td_var%t_dim(2)%i_len, & 
     596         &                td_var%t_dim(3)%i_len, & 
     597         &                td_var%t_dim(4)%i_len) ) 
     598 
     599         bl_mask(:,:,:,:)=1 
     600         WHERE(td_var%d_value(:,:,:,:)==td_var%d_fill) bl_mask(:,:,:,:)=0       
     601 
     602         SELECT CASE(TRIM(td_var%c_point)) 
     603         CASE DEFAULT ! 'T' 
     604            tl_mask=var_init('tmask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 
     605            &                id_ew=td_var%i_ew ) 
     606         CASE('U') 
     607            tl_mask=var_init('umask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 
     608            &                id_ew=td_var%i_ew ) 
     609         CASE('V') 
     610            tl_mask=var_init('vmask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 
     611            &                id_ew=td_var%i_ew ) 
     612         CASE('F') 
     613            tl_mask=var_init('fmask',bl_mask(:,:,:,:),td_dim=td_var%t_dim(:),& 
     614            &                id_ew=td_var%i_ew ) 
     615         END SELECT          
     616 
     617         DEALLOCATE(bl_mask) 
     618 
     619         ! interpolate mask 
     620         CALL interp_fill_value( tl_mask, id_rho(:), & 
     621         &                       id_offset=id_offset(:,:) ) 
     622 
     623         ! work on variable 
     624         ! add extraband 
     625         CALL extrap_add_extrabands(td_var, il_iext, il_jext) 
     626 
     627         ! extrapolate variable 
     628         CALL extrap_fill_value( td_var, id_iext=il_iext, id_jext=il_jext ) 
     629 
     630         ! interpolate variable 
     631         CALL interp_fill_value( td_var, id_rho(:), & 
     632         &                       id_offset=id_offset(:,:)) 
     633 
     634         ! remove extraband 
     635         CALL extrap_del_extrabands(td_var, il_iext*id_rho(jp_I), il_jext*id_rho(jp_J)) 
     636 
     637         ! keep original mask  
     638         WHERE( tl_mask%d_value(:,:,:,:) == 0 ) 
     639            td_var%d_value(:,:,:,:)=td_var%d_fill 
     640         END WHERE 
     641      ENDIF 
    443642 
    444643      ! clean variable structure 
    445644      CALL var_clean(tl_mask) 
    446       CALL var_clean(tl_var) 
    447645 
    448646   END SUBROUTINE create_coord_interp 
    449    !> @endcode 
    450647END PROGRAM create_coord 
Note: See TracChangeset for help on using the changeset viewer.