Ignore:
Timestamp:
2020-02-25T16:29:34+01:00 (13 months ago)
Author:
jcastill
Message:

First implementation of the branch - compiling after merge

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/DOM/dtatsd.F90

    r11715 r12453  
    2121   ! 
    2222   USE in_out_manager  ! I/O manager 
     23   USE iom             ! IOM library 
    2324   USE lib_mpp         ! MPP library 
    2425 
     
    3132   !                                  !!* namtsd  namelist : Temperature & Salinity Data * 
    3233   LOGICAL , PUBLIC ::   ln_tsd_init   !: T & S data flag 
     34   LOGICAL , PUBLIC ::   ln_tsd_interp !: vertical interpolation flag 
    3335   LOGICAL , PUBLIC ::   ln_tsd_dmp    !: internal damping toward input data flag 
    3436 
    3537   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_tsd   ! structure of input SST (file informations, fields read) 
     38   INTEGER                              ::   jpk_init, inum_dta  
     39   INTEGER                              ::   id, linum   ! local integers  
     40   INTEGER                              ::   zdim(4) 
    3641 
    3742   !!---------------------------------------------------------------------- 
     
    5358      LOGICAL, INTENT(in), OPTIONAL ::   ld_tradmp   ! force the initialization when tradp is used 
    5459      ! 
    55       INTEGER ::   ios, ierr0, ierr1, ierr2, ierr3   ! local integers 
     60      INTEGER ::   ios, ierr0, ierr1, ierr2, ierr3, ierr4, ierr5   ! local integers 
    5661      !! 
    5762      CHARACTER(len=100)            ::   cn_dir          ! Root directory for location of ssr files 
    58       TYPE(FLD_N), DIMENSION( jpts) ::   slf_i           ! array of namelist informations on the fields to read 
    59       TYPE(FLD_N)                   ::   sn_tem, sn_sal 
     63      TYPE(FLD_N), DIMENSION(jpts+2)::   slf_i           ! array of namelist informations on the fields to read  
     64      TYPE(FLD_N)                   ::   sn_tem, sn_sal, sn_dep, sn_msk 
    6065      !! 
    61       NAMELIST/namtsd/   ln_tsd_init, ln_tsd_dmp, cn_dir, sn_tem, sn_sal 
     66      NAMELIST/namtsd/   ln_tsd_init, ln_tsd_interp, ln_tsd_dmp, cn_dir, sn_tem, sn_sal, sn_dep, sn_msk 
    6267      !!---------------------------------------------------------------------- 
    6368      ! 
    6469      !  Initialisation 
    65       ierr0 = 0  ;  ierr1 = 0  ;  ierr2 = 0  ;  ierr3 = 0 
     70      ierr0 = 0  ;  ierr1 = 0  ;  ierr2 = 0  ;  ierr3 = 0  ; ierr4 = 0  ;  ierr5 = 0 
    6671      ! 
    6772      REWIND( numnam_ref )              ! Namelist namtsd in reference namelist :  
     
    8085         WRITE(numout,*) '~~~~~~~~~~~~ ' 
    8186         WRITE(numout,*) '   Namelist namtsd' 
    82          WRITE(numout,*) '      Initialisation of ocean T & S with T &S input data   ln_tsd_init = ', ln_tsd_init 
    83          WRITE(numout,*) '      damping of ocean T & S toward T &S input data        ln_tsd_dmp  = ', ln_tsd_dmp 
     87         WRITE(numout,*) '      Initialisation of ocean T & S with T &S input data   ln_tsd_init   = ', ln_tsd_init 
     88         WRITE(numout,*) '      Interpolation of initial conditions in the vertical  ln_tsd_interp = ', ln_tsd_interp 
     89         WRITE(numout,*) '      damping of ocean T & S toward T &S input data        ln_tsd_dmp    = ', ln_tsd_dmp 
    8490         WRITE(numout,*) 
    8591         IF( .NOT.ln_tsd_init .AND. .NOT.ln_tsd_dmp ) THEN 
     
    94100         ln_tsd_init = .FALSE. 
    95101      ENDIF 
     102      IF( ln_tsd_interp .AND. ln_tsd_dmp ) THEN  
     103            CALL ctl_stop( 'dta_tsd_init: Tracer damping and vertical interpolation not yet configured' )   ;   RETURN  
     104      ENDIF  
     105      IF( ln_tsd_interp .AND. LEN(TRIM(sn_msk%wname)) > 0 ) THEN  
     106            CALL ctl_stop( 'dta_tsd_init: Using vertical interpolation and weights files not recommended' )   ;   RETURN  
     107      ENDIF 
    96108      ! 
    97109      !                             ! allocate the arrays (if necessary) 
    98110      IF( ln_tsd_init .OR. ln_tsd_dmp ) THEN 
    99111         ! 
    100          ALLOCATE( sf_tsd(jpts), STAT=ierr0 ) 
     112         IF( ln_tsd_interp ) THEN  
     113           ALLOCATE( sf_tsd(jpts+2), STAT=ierr0 ) ! to carry the addtional depth information  
     114         ELSE  
     115           ALLOCATE( sf_tsd(jpts  ), STAT=ierr0 )   
     116         ENDIF 
    101117         IF( ierr0 > 0 ) THEN 
    102118            CALL ctl_stop( 'dta_tsd_init: unable to allocate sf_tsd structure' )   ;   RETURN 
    103119         ENDIF 
    104120         ! 
    105                                 ALLOCATE( sf_tsd(jp_tem)%fnow(jpi,jpj,jpk)   , STAT=ierr0 ) 
    106          IF( sn_tem%ln_tint )   ALLOCATE( sf_tsd(jp_tem)%fdta(jpi,jpj,jpk,2) , STAT=ierr1 ) 
    107                                 ALLOCATE( sf_tsd(jp_sal)%fnow(jpi,jpj,jpk)   , STAT=ierr2 ) 
    108          IF( sn_sal%ln_tint )   ALLOCATE( sf_tsd(jp_sal)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 ) 
    109          ! 
    110          IF( ierr0 + ierr1 + ierr2 + ierr3 > 0 ) THEN 
     121         IF( ln_tsd_interp ) THEN  
     122            CALL iom_open ( trim(cn_dir) // trim(sn_dep%clname), inum_dta )   
     123            id = iom_varid( inum_dta, sn_dep%clvar, zdim )  
     124            jpk_init = zdim(3)  
     125            IF(lwp) WRITE(numout,*) 'Dimension of vertical coordinate in ICs: ', jpk_init  
     126            CALL iom_close( inum_dta )   ! Close the input file  
     127            !  
     128                                 ALLOCATE( sf_tsd(jp_tem)%fnow(jpi,jpj,jpk_init  ) , STAT=ierr0 )  
     129            IF( sn_tem%ln_tint ) ALLOCATE( sf_tsd(jp_tem)%fdta(jpi,jpj,jpk_init,2) , STAT=ierr1 )  
     130                                 ALLOCATE( sf_tsd(jp_sal)%fnow(jpi,jpj,jpk_init  ) , STAT=ierr2 )  
     131            IF( sn_sal%ln_tint ) ALLOCATE( sf_tsd(jp_sal)%fdta(jpi,jpj,jpk_init,2) , STAT=ierr3 )    
     132                                 ALLOCATE( sf_tsd(jp_dep)%fnow(jpi,jpj,jpk_init  ) , STAT=ierr4 )  
     133                                 ALLOCATE( sf_tsd(jp_msk)%fnow(jpi,jpj,jpk_init  ) , STAT=ierr5 )  
     134         ELSE  
     135                                 ALLOCATE( sf_tsd(jp_tem)%fnow(jpi,jpj,jpk)   , STAT=ierr0 )  
     136            IF( sn_tem%ln_tint ) ALLOCATE( sf_tsd(jp_tem)%fdta(jpi,jpj,jpk,2) , STAT=ierr1 )  
     137                                 ALLOCATE( sf_tsd(jp_sal)%fnow(jpi,jpj,jpk)   , STAT=ierr2 )  
     138            IF( sn_sal%ln_tint ) ALLOCATE( sf_tsd(jp_sal)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 )    
     139         ENDIF ! ln_tsd_interp  
     140  
     141         !  
     142         IF( ierr0 + ierr1 + ierr2 + ierr3 + ierr4 + ierr5 > 0 ) THEN 
    111143            CALL ctl_stop( 'dta_tsd : unable to allocate T & S data arrays' )   ;   RETURN 
    112144         ENDIF 
    113145         !                         ! fill sf_tsd with sn_tem & sn_sal and control print 
    114146         slf_i(jp_tem) = sn_tem   ;   slf_i(jp_sal) = sn_sal 
     147         IF( ln_tsd_interp ) slf_i(jp_dep) = sn_dep   ;   slf_i(jp_msk) = sn_msk 
    115148         CALL fld_fill( sf_tsd, slf_i, cn_dir, 'dta_tsd', 'Temperature & Salinity data', 'namtsd', no_print ) 
    116149         ! 
     
    138171      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(  out) ::   ptsd   ! T & S data 
    139172      ! 
    140       INTEGER ::   ji, jj, jk, jl, jkk   ! dummy loop indicies 
    141       INTEGER ::   ik, il0, il1, ii0, ii1, ij0, ij1   ! local integers 
     173      INTEGER ::   ji, jj, jk, jl, jk_init   ! dummy loop indices  
     174      INTEGER ::   ik, il0, il1, ii0, ii1, ij0, ij1        ! local integers 
    142175      REAL(wp)::   zl, zi                             ! local scalars 
    143       REAL(wp), DIMENSION(jpk) ::  ztp, zsp   ! 1D workspace 
    144176      !!---------------------------------------------------------------------- 
    145177      ! 
     
    176208!!gm end 
    177209      ! 
    178       ptsd(:,:,:,jp_tem) = sf_tsd(jp_tem)%fnow(:,:,:)    ! NO mask 
    179       ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:)  
    180       ! 
    181       IF( ln_sco ) THEN                   !==   s- or mixed s-zps-coordinate   ==! 
    182          ! 
    183          IF( kt == nit000 .AND. lwp )THEN 
    184             WRITE(numout,*) 
    185             WRITE(numout,*) 'dta_tsd: interpolates T & S data onto the s- or mixed s-z-coordinate mesh' 
    186          ENDIF 
    187          ! 
    188          DO jj = 1, jpj                         ! vertical interpolation of T & S 
    189             DO ji = 1, jpi 
    190                DO jk = 1, jpk                        ! determines the intepolated T-S profiles at each (i,j) points 
     210      IF( kt == nit000 .AND. lwp )THEN  
     211         WRITE(numout,*)  
     212         WRITE(numout,*) 'dta_tsd: interpolates T & S data onto current mesh'  
     213      ENDIF  
     214      !  
     215      IF( ln_tsd_interp ) THEN                 ! probably should use pointers in the following to make more readable  
     216      !  
     217         DO jk = 1, jpk                        ! determines the intepolated T-S profiles at each (i,j) points  
     218            DO jj= 1, jpj  
     219               DO ji= 1, jpi 
    191220                  zl = gdept_0(ji,jj,jk) 
    192                   IF(     zl < gdept_1d(1  ) ) THEN          ! above the first level of data 
    193                      ztp(jk) =  ptsd(ji,jj,1    ,jp_tem) 
    194                      zsp(jk) =  ptsd(ji,jj,1    ,jp_sal) 
    195                   ELSEIF( zl > gdept_1d(jpk) ) THEN          ! below the last level of data 
    196                      ztp(jk) =  ptsd(ji,jj,jpkm1,jp_tem) 
    197                      zsp(jk) =  ptsd(ji,jj,jpkm1,jp_sal) 
    198                   ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1 
    199                      DO jkk = 1, jpkm1                                  ! when  gdept(jkk) < zl < gdept(jkk+1) 
    200                         IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN 
    201                            zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk)) 
    202                            ztp(jk) = ptsd(ji,jj,jkk,jp_tem) + ( ptsd(ji,jj,jkk+1,jp_tem) - ptsd(ji,jj,jkk,jp_tem) ) * zi  
    203                            zsp(jk) = ptsd(ji,jj,jkk,jp_sal) + ( ptsd(ji,jj,jkk+1,jp_sal) - ptsd(ji,jj,jkk,jp_sal) ) * zi 
     221                  IF( zl < sf_tsd(jp_dep)%fnow(ji,jj,1) ) THEN                     ! above the first level of data  
     222                     ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,1)   
     223                     ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,1)  
     224                  ELSEIF( zl > sf_tsd(jp_dep)%fnow(ji,jj,jpk_init) ) THEN          ! below the last level of data  
     225                     ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,jpk_init)  
     226                     ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,jpk_init)  
     227                  ELSE                                                             ! inbetween : vertical interpolation between jk_init & jk_init+1  
     228                     DO jk_init = 1, jpk_init-1                                    ! when  gdept(jk_init) < zl < gdept(jk_init+1)  
     229                        IF( sf_tsd(jp_msk)%fnow(ji,jj,jk_init+1) == 0 ) THEN       ! if there is no data fill down  
     230                           sf_tsd(jp_tem)%fnow(ji,jj,jk_init+1) = sf_tsd(jp_tem)%fnow(ji,jj,jk_init)  
     231                           sf_tsd(jp_sal)%fnow(ji,jj,jk_init+1) = sf_tsd(jp_sal)%fnow(ji,jj,jk_init)  
     232                        ENDIF  
     233                        IF((zl-sf_tsd(jp_dep)%fnow(ji,jj,jk_init)) * (zl-sf_tsd(jp_dep)%fnow(ji,jj,jk_init+1)) <= 0._wp ) THEN  
     234                           zi = ( zl - sf_tsd(jp_dep)%fnow(ji,jj,jk_init) ) / &  
     235                        &       (sf_tsd(jp_dep)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_dep)%fnow(ji,jj,jk_init))  
     236                           ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,jk_init) + &  
     237                        &                          (sf_tsd(jp_tem)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_tem)%fnow(ji,jj,jk_init)) * zi  
     238                           ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,jk_init) + &  
     239                        &                          (sf_tsd(jp_sal)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_sal)%fnow(ji,jj,jk_init)) * zi 
    204240                        ENDIF 
    205241                     END DO 
    206242                  ENDIF 
    207243               END DO 
    208                DO jk = 1, jpkm1 
    209                   ptsd(ji,jj,jk,jp_tem) = ztp(jk) * tmask(ji,jj,jk)     ! mask required for mixed zps-s-coord 
    210                   ptsd(ji,jj,jk,jp_sal) = zsp(jk) * tmask(ji,jj,jk) 
    211                END DO 
    212                ptsd(ji,jj,jpk,jp_tem) = 0._wp 
    213                ptsd(ji,jj,jpk,jp_sal) = 0._wp 
    214244            END DO 
    215245         END DO 
    216246         !  
     247         ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) *tmask(:,:,:)  
     248         ptsd(:,:,:,jp_sal) = ptsd(:,:,:,jp_sal) *tmask(:,:,:) 
    217249      ELSE                                !==   z- or zps- coordinate   ==! 
    218250         !                              
    219          ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) * tmask(:,:,:)    ! Mask 
    220          ptsd(:,:,:,jp_sal) = ptsd(:,:,:,jp_sal) * tmask(:,:,:) 
     251         ptsd(:,:,:,jp_tem) = sf_tsd(jp_tem)%fnow(:,:,:)  * tmask(:,:,:)  ! Mask  
     252         ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:) * tmask(:,:,:) 
    221253         ! 
    222254         IF( ln_zps ) THEN                      ! zps-coordinate (partial steps) interpolation at the last ocean level 
     
    248280                                        DEALLOCATE( sf_tsd(jp_sal)%fnow )     ! S arrays in the structure 
    249281         IF( sf_tsd(jp_sal)%ln_tint )   DEALLOCATE( sf_tsd(jp_sal)%fdta ) 
     282         IF( ln_tsd_interp )            DEALLOCATE( sf_tsd(jp_dep)%fnow )     ! T arrays in the structure  
     283         IF( ln_tsd_interp )            DEALLOCATE( sf_tsd(jp_msk)%fnow )     ! T arrays in the structure  
    250284                                        DEALLOCATE( sf_tsd              )     ! the structure itself 
    251285      ENDIF 
Note: See TracChangeset for help on using the changeset viewer.