New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 7773 for branches/UKMO/dev_r5785_SSS_obsoper/NEMOGCM/NEMO/OPA_SRC/OBS/obs_sstbias.F90 – NEMO

Ignore:
Timestamp:
2017-03-09T13:52:43+01:00 (7 years ago)
Author:
mattmartin
Message:

Committing updates after doing the following:

  • merging the branch dev_r4650_general_vert_coord_obsoper@7763 into this branch
  • updating it so that the following OBS changes were implemented correctly on top of the simplification changes:
    • generalised vertical coordinate for profile obs. This was done so that is now the default option.
    • sst bias correction implemented with the new simplified obs code.
    • included the biogeochemical obs types int he new simplified obs code.
    • included the changes to exclude obs in the boundary for limited area models
    • included other changes for the efficiency of the obs operator to remove global arrays.
File:
1 copied

Legend:

Unmodified
Added
Removed

  • branches/UKMO/dev_r5785_SSS_obsoper/NEMOGCM/NEMO/OPA_SRC/OBS/obs_sstbias.F90

    r7740 r7773  
    11MODULE obs_sstbias 
    22   !!====================================================================== 
    3    !!                       ***  MODULE obs_readsstbias  *** 
    4    !! Observation diagnostics: Read the bias for SLA data 
     3   !!                       ***  MODULE obs_sstbias  *** 
     4   !! Observation diagnostics: Read the bias for SST data 
    55   !!====================================================================== 
    66   !!---------------------------------------------------------------------- 
    7    !!   obs_rea_sstbias : Driver for reading altimeter bias 
     7   !!   obs_app_sstbias : Driver for reading and applying the SST bias 
    88   !!---------------------------------------------------------------------- 
    99   !! * Modules used    
    … …  
    2222   USE dom_oce, ONLY : &        ! Domain variables 
    2323      & tmask, & 
    24       & tmask_i, & 
    25       & e1t,   & 
    26       & e2t,   & 
    2724      & gphit, & 
    2825      & glamt 
    29    USE oce, ONLY : &           ! Model variables 
    30       & sshn 
    3126   USE obs_inter_h2d 
    3227   USE obs_utils               ! Various observation tools 
    … …  
    3732   PUBLIC obs_app_sstbias     ! Read the altimeter bias 
    3833CONTAINS 
    39    SUBROUTINE obs_app_sstbias( ksstno, sstdata, k2dint, knumtypes, & 
     34   SUBROUTINE obs_app_sstbias( sstdata, k2dint, knumtypes, & 
    4035                               cl_bias_files ) 
    4136      !!--------------------------------------------------------------------- 
    4237      !! 
    43       !!                   *** ROUTINE obs_rea_sstbias *** 
     38      !!                   *** ROUTINE obs_app_sstbias *** 
    4439      !! 
    4540      !! ** Purpose : Read SST bias data from files and apply correction to 
    … …  
    5954      USE iom 
    6055      USE netcdf 
     56 
    6157      !! * Arguments 
    62       INTEGER, INTENT(IN) :: ksstno      ! Number of SST obs sets 
    63       TYPE(obs_surf), DIMENSION(ksstno), INTENT(INOUT) :: & 
    64          & sstdata       ! SST data 
     58      TYPE(obs_surf), INTENT(INOUT) :: & 
     59         & sstdata            ! SST data 
    6560      INTEGER, INTENT(IN) :: k2dint 
    66       INTEGER, INTENT(IN) :: knumtypes !number of bias types to read in 
     61      INTEGER, INTENT(IN) :: & 
     62         & knumtypes          ! Number of bias types to read in 
    6763      CHARACTER(LEN=128), DIMENSION(knumtypes), INTENT(IN) :: & 
    68                           cl_bias_files !List of files to read 
     64         & cl_bias_files      ! List of files to read 
     65 
    6966      !! * Local declarations 
    7067      INTEGER :: jslano       ! Data set loop variable 
    … …  
    8077      INTEGER :: i_var_id 
    8178      INTEGER, DIMENSION(knumtypes) :: & 
    82          & ibiastypes             ! Array of the bias types in each file 
     79         & ibiastypes         ! Array of the bias types in each file 
    8380      REAL(wp), DIMENSION(jpi,jpj,knumtypes) :: &  
    84          & z_sstbias              ! Array to store the SST bias values 
     81         & z_sstbias          ! Array to store the SST bias values 
    8582      REAL(wp), DIMENSION(jpi,jpj) :: &  
    86          & z_sstbias_2d           ! Array to store the SST bias values    
     83         & z_sstbias_2d       ! Array to store the SST bias values    
    8784      REAL(wp), DIMENSION(1) :: & 
    8885         & zext, & 
    … …  
    114111      INTEGER :: iret  
    115112      INTEGER :: inumtype 
    116       IF(lwp)WRITE(numout,*)  
    117       IF(lwp)WRITE(numout,*) 'obs_rea_sstbias : ' 
    118       IF(lwp)WRITE(numout,*) '----------------- ' 
    119       IF(lwp)WRITE(numout,*) 'Read SST bias ' 
    120       ! Open and read the files 
    121       z_sstbias(:,:,:)=0.0_wp 
     113 
     114      IF ( lwp ) THEN 
     115         WRITE(numout,*)  
     116         WRITE(numout,*) 'obs_app_sstbias : ' 
     117         WRITE(numout,*) '----------------- ' 
     118         WRITE(numout,*) 'Read SST bias ' 
     119      ENDIF 
     120 
     121      ! Open and read the SST bias files for each bias type 
     122      z_sstbias(:,:,:) = 0.0_wp 
     123 
    122124      DO jtype = 1, knumtypes 
    123125      
    124126         numsstbias=0 
    125          IF(lwp)WRITE(numout,*) 'Opening ',cl_bias_files(jtype) 
    126          CALL iom_open( cl_bias_files(jtype), numsstbias, ldstop=.FALSE. )        
     127 
     128         IF ( lwp ) WRITE(numout,*) 'Opening ',cl_bias_files(jtype) 
     129         CALL iom_open( cl_bias_files(jtype), numsstbias, ldstop=.FALSE. ) 
     130 
    127131         IF (numsstbias .GT. 0) THEN 
    128132      
    … …  
    137141            iret=NF90_CLOSE(incfile)        
    138142            
    139             IF ( iret /= 0  ) CALL ctl_stop( & 
    140                'obs_rea_sstbias : Cannot read bias type from file '// & 
    141                cl_bias_files(jtype) ) 
     143            IF ( iret /= 0  ) THEN 
     144               CALL ctl_stop( 'obs_app_sstbias : Cannot read bias type from file '// & 
     145                  &           TRIM( cl_bias_files(jtype) ) ) 
     146            ENDIF 
     147 
    142148            ! Get the SST bias data 
    143149            CALL iom_get( numsstbias, jpdom_data, 'tn', z_sstbias_2d(:,:), 1 ) 
    144150            z_sstbias(:,:,jtype) = z_sstbias_2d(:,:)        
    145151            ! Close the file 
    146             CALL iom_close(numsstbias)        
     152            CALL iom_close(numsstbias) 
     153      
    147154         ELSE 
    148155            CALL ctl_stop('obs_read_sstbias: File '// &  
    149                            TRIM( cl_bias_files(jtype) )//' Not found') 
     156               &          TRIM( cl_bias_files(jtype) )//' Not found') 
    150157         ENDIF 
     158 
    151159      END DO 
    152160            
    153       ! Interpolate the bias already on the model grid at the observation point 
    154       DO jslano = 1, ksstno 
     161      ! Interpolate the bias from the model grid to the observation points 
     162      ALLOCATE( & 
     163         & igrdi(2,2,sstdata%nsurf), & 
     164         & igrdj(2,2,sstdata%nsurf), & 
     165         & zglam(2,2,sstdata%nsurf), & 
     166         & zgphi(2,2,sstdata%nsurf), & 
     167         & zmask(2,2,sstdata%nsurf)  ) 
     168        
     169      DO jobs = 1, sstdata%nsurf  
     170         igrdi(1,1,jobs) = sstdata%mi(jobs)-1 
     171         igrdj(1,1,jobs) = sstdata%mj(jobs)-1 
     172         igrdi(1,2,jobs) = sstdata%mi(jobs)-1 
     173         igrdj(1,2,jobs) = sstdata%mj(jobs) 
     174         igrdi(2,1,jobs) = sstdata%mi(jobs) 
     175         igrdj(2,1,jobs) = sstdata%mj(jobs)-1 
     176         igrdi(2,2,jobs) = sstdata%mi(jobs) 
     177         igrdj(2,2,jobs) = sstdata%mj(jobs) 
     178      END DO 
     179 
     180      CALL obs_int_comm_2d( 2, 2, sstdata%nsurf, & 
     181         &                  igrdi, igrdj, glamt, zglam ) 
     182      CALL obs_int_comm_2d( 2, 2, sstdata%nsurf, & 
     183         &                  igrdi, igrdj, gphit, zgphi ) 
     184      CALL obs_int_comm_2d( 2, 2, sstdata%nsurf, & 
     185         &                  igrdi, igrdj, tmask(:,:,1), zmask ) 
     186 
     187      DO jtype = 1, knumtypes 
     188          
     189         !Find the number observations of type 
     190         !and alllocate tempory arrays 
     191         inumtype = COUNT( sstdata%ntyp(:) == ibiastypes(jtype) ) 
     192 
    155193         ALLOCATE( & 
    156             & igrdi(2,2,sstdata(jslano)%nsurf), & 
    157             & igrdj(2,2,sstdata(jslano)%nsurf), & 
    158             & zglam(2,2,sstdata(jslano)%nsurf), & 
    159             & zgphi(2,2,sstdata(jslano)%nsurf), & 
    160             & zmask(2,2,sstdata(jslano)%nsurf)  ) 
    161         
    162          DO jobs = 1, sstdata(jslano)%nsurf  
    163             igrdi(1,1,jobs) = sstdata(jslano)%mi(jobs)-1 
    164             igrdj(1,1,jobs) = sstdata(jslano)%mj(jobs)-1 
    165             igrdi(1,2,jobs) = sstdata(jslano)%mi(jobs)-1 
    166             igrdj(1,2,jobs) = sstdata(jslano)%mj(jobs) 
    167             igrdi(2,1,jobs) = sstdata(jslano)%mi(jobs) 
    168             igrdj(2,1,jobs) = sstdata(jslano)%mj(jobs)-1 
    169             igrdi(2,2,jobs) = sstdata(jslano)%mi(jobs) 
    170             igrdj(2,2,jobs) = sstdata(jslano)%mj(jobs) 
    171          END DO 
    172          CALL obs_int_comm_2d( 2, 2, sstdata(jslano)%nsurf, & 
    173             &                  igrdi, igrdj, glamt, zglam ) 
    174          CALL obs_int_comm_2d( 2, 2, sstdata(jslano)%nsurf, & 
    175             &                  igrdi, igrdj, gphit, zgphi ) 
    176          CALL obs_int_comm_2d( 2, 2, sstdata(jslano)%nsurf, & 
    177             &                  igrdi, igrdj, tmask(:,:,1), zmask ) 
    178          DO jtype = 1, knumtypes 
    179           
    180             !Find the number observations of type 
    181             !and alllocate tempory arrays 
    182             inumtype = COUNT( sstdata(jslano)%ntyp(:) == ibiastypes(jtype) ) 
    183             ALLOCATE( & 
    184194            & igrdi_tmp(2,2,inumtype), & 
    185195            & igrdj_tmp(2,2,inumtype), & 
    … …  
    188198            & zmask_tmp(2,2,inumtype), & 
    189199            & zbias( 2,2,inumtype ) ) 
    190             jt=1 
    191             DO jobs = 1, sstdata(jslano)%nsurf  
    192                IF ( sstdata(jslano)%ntyp(jobs) == ibiastypes(jtype) ) THEN 
    193                   igrdi_tmp(:,:,jt) = igrdi(:,:,jobs)  
    194                   igrdj_tmp(:,:,jt) = igrdj(:,:,jobs) 
    195                   zglam_tmp(:,:,jt) = zglam(:,:,jobs) 
    196                   zgphi_tmp(:,:,jt) = zgphi(:,:,jobs) 
    197                   zgphi_tmp(:,:,jt) = zgphi(:,:,jobs) 
    198                   zmask_tmp(:,:,jt) = zmask(:,:,jobs) 
    199                   jt = jt +1 
    200                ENDIF 
    201             END DO 
     200 
     201         jt=1 
     202         DO jobs = 1, sstdata%nsurf  
     203            IF ( sstdata%ntyp(jobs) == ibiastypes(jtype) ) THEN 
     204 
     205               igrdi_tmp(:,:,jt) = igrdi(:,:,jobs)  
     206               igrdj_tmp(:,:,jt) = igrdj(:,:,jobs) 
     207               zglam_tmp(:,:,jt) = zglam(:,:,jobs) 
     208               zgphi_tmp(:,:,jt) = zgphi(:,:,jobs) 
     209               zgphi_tmp(:,:,jt) = zgphi(:,:,jobs) 
     210               zmask_tmp(:,:,jt) = zmask(:,:,jobs) 
     211 
     212               jt = jt +1 
     213 
     214            ENDIF 
     215         END DO 
    202216                          
    203             CALL obs_int_comm_2d( 2, 2, inumtype, & 
    204                   &           igrdi_tmp(:,:,:), igrdj_tmp(:,:,:), & 
    205                   &           z_sstbias(:,:,jtype), zbias(:,:,:) ) 
    206             jt=1 
    207             DO jobs = 1, sstdata(jslano)%nsurf 
    208                IF ( sstdata(jslano)%ntyp(jobs) == ibiastypes(jtype) ) THEN 
    209                   zlam = sstdata(jslano)%rlam(jobs) 
    210                   zphi = sstdata(jslano)%rphi(jobs) 
    211                   iico = sstdata(jslano)%mi(jobs) 
    212                   ijco = sstdata(jslano)%mj(jobs)          
    213                   CALL obs_int_h2d_init( 1, 1, k2dint, zlam, zphi,         & 
    214                      &                   zglam_tmp(:,:,jt), & 
    215                      &                   zgphi_tmp(:,:,jt), & 
    216                      &                   zmask_tmp(:,:,jt), zweig, zobsmask ) 
    217                   CALL obs_int_h2d( 1, 1,      & 
    218                      &              zweig, zbias(:,:,jt),  zext ) 
    219                   ! adjust sst with bias field 
    220                   sstdata(jslano)%robs(jobs,1) = & 
    221                      sstdata(jslano)%robs(jobs,1) - zext(1) 
    222                   jt=jt+1 
    223                ENDIF 
    224             END DO  
     217         CALL obs_int_comm_2d( 2, 2, inumtype, & 
     218               &           igrdi_tmp(:,:,:), igrdj_tmp(:,:,:), & 
     219               &           z_sstbias(:,:,jtype), zbias(:,:,:) ) 
     220 
     221         jt=1 
     222         DO jobs = 1, sstdata%nsurf 
     223            IF ( sstdata%ntyp(jobs) == ibiastypes(jtype) ) THEN 
     224 
     225               zlam = sstdata%rlam(jobs) 
     226               zphi = sstdata%rphi(jobs) 
     227               iico = sstdata%mi(jobs) 
     228               ijco = sstdata%mj(jobs)          
     229 
     230               CALL obs_int_h2d_init( 1, 1, k2dint, zlam, zphi,         & 
     231                  &                   zglam_tmp(:,:,jt), & 
     232                  &                   zgphi_tmp(:,:,jt), & 
     233                  &                   zmask_tmp(:,:,jt), zweig, zobsmask ) 
     234 
     235               CALL obs_int_h2d( 1, 1, zweig, zbias(:,:,jt),  zext ) 
     236 
     237               ! adjust sst with bias field 
     238               sstdata%robs(jobs,1) = & 
     239                  &    sstdata%robs(jobs,1) - zext(1) 
     240 
     241               jt=jt+1 
     242 
     243            ENDIF 
     244         END DO  
    225245                
    226             !Deallocate arrays 
    227             DEALLOCATE( & 
    228             & igrdi_tmp, & 
    229             & igrdj_tmp, & 
    230             & zglam_tmp, & 
    231             & zgphi_tmp, & 
    232             & zmask_tmp, & 
    233             & zbias )            
    234          END DO 
     246         !Deallocate arrays 
    235247         DEALLOCATE( & 
    236             & igrdi, & 
    237             & igrdj, & 
    238             & zglam, & 
    239             & zgphi, & 
    240             & zmask ) 
    241       END DO 
     248         & igrdi_tmp, & 
     249         & igrdj_tmp, & 
     250         & zglam_tmp, & 
     251         & zgphi_tmp, & 
     252         & zmask_tmp, & 
     253         & zbias )       
     254      
     255      END DO !jtype 
     256 
     257      DEALLOCATE( & 
     258         & igrdi, & 
     259         & igrdj, & 
     260         & zglam, & 
     261         & zgphi, & 
     262         & zmask ) 
     263 
    242264      IF(lwp) THEN 
    243265         WRITE(numout,*) " " 
    244266         WRITE(numout,*) "SST bias correction applied successfully" 
    245267         WRITE(numout,*) "Obs types: ",ibiastypes(:), & 
    246                               " Have all been bias corrected\n" 
     268                              " have all been bias corrected\n" 
    247269      ENDIF 
    248270   END SUBROUTINE obs_app_sstbias 
Note: See TracChangeset for help on using the changeset viewer.