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 9023 for branches/2017/dev_merge_2017/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90 – NEMO

Ignore:
Timestamp:
2017-12-13T18:08:50+01:00 (6 years ago)
Author:
timgraham
Message:

Merged METO_MERCATOR branch and resolved all conflicts in OPA_SRC

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_merge_2017/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90

    r9019 r9023  
    4646   LOGICAL, PUBLIC :: ln_diaobs   !: Logical switch for the obs operator 
    4747   LOGICAL :: ln_sstnight         !: Logical switch for night mean SST obs 
    48     
     48   LOGICAL :: ln_sla_fp_indegs    !: T=> SLA obs footprint size specified in degrees, F=> in metres 
     49   LOGICAL :: ln_sst_fp_indegs    !: T=> SST obs footprint size specified in degrees, F=> in metres 
     50   LOGICAL :: ln_sss_fp_indegs    !: T=> SSS obs footprint size specified in degrees, F=> in metres 
     51   LOGICAL :: ln_sic_fp_indegs    !: T=> sea-ice obs footprint size specified in degrees, F=> in metres 
     52 
     53   REAL(wp) :: rn_sla_avglamscl !: E/W diameter of SLA observation footprint (metres) 
     54   REAL(wp) :: rn_sla_avgphiscl !: N/S diameter of SLA observation footprint (metres) 
     55   REAL(wp) :: rn_sst_avglamscl !: E/W diameter of SST observation footprint (metres) 
     56   REAL(wp) :: rn_sst_avgphiscl !: N/S diameter of SST observation footprint (metres) 
     57   REAL(wp) :: rn_sss_avglamscl !: E/W diameter of SSS observation footprint (metres) 
     58   REAL(wp) :: rn_sss_avgphiscl !: N/S diameter of SSS observation footprint (metres) 
     59   REAL(wp) :: rn_sic_avglamscl !: E/W diameter of sea-ice observation footprint (metres) 
     60   REAL(wp) :: rn_sic_avgphiscl !: N/S diameter of sea-ice observation footprint (metres) 
     61 
    4962   INTEGER :: nn_1dint       !: Vertical interpolation method 
    50    INTEGER :: nn_2dint       !: Horizontal interpolation method 
     63   INTEGER :: nn_2dint       !: Default horizontal interpolation method 
     64   INTEGER :: nn_2dint_sla   !: SLA horizontal interpolation method  
     65   INTEGER :: nn_2dint_sst   !: SST horizontal interpolation method  
     66   INTEGER :: nn_2dint_sss   !: SSS horizontal interpolation method  
     67   INTEGER :: nn_2dint_sic   !: Seaice horizontal interpolation method  
    5168   INTEGER, DIMENSION(imaxavtypes) ::   nn_profdavtypes   !: Profile data types representing a daily average 
    5269   INTEGER :: nproftypes     !: Number of profile obs types 
    5370   INTEGER :: nsurftypes     !: Number of surface obs types 
    54    INTEGER, DIMENSION(:), ALLOCATABLE ::   nvarsprof, nvarssurf   !: Number of profile & surface variables 
    55    INTEGER, DIMENSION(:), ALLOCATABLE ::   nextrprof, nextrsurf   !: Number of profile & surface extra variables 
    56    INTEGER, PUBLIC, ALLOCATABLE, DIMENSION(:) ::   sstbias_type   !: SST bias type     
    57    TYPE(obs_surf), PUBLIC, POINTER, DIMENSION(:) ::   surfdata, surfdataqc   !: Initial surface data before & after quality control 
    58    TYPE(obs_prof), PUBLIC, POINTER, DIMENSION(:) ::   profdata, profdataqc   !: Initial profile data before & after quality control 
     71   INTEGER, DIMENSION(:), ALLOCATABLE :: & 
     72      & nvarsprof, &         !: Number of profile variables 
     73      & nvarssurf            !: Number of surface variables 
     74   INTEGER, DIMENSION(:), ALLOCATABLE :: & 
     75      & nextrprof, &         !: Number of profile extra variables 
     76      & nextrsurf            !: Number of surface extra variables 
     77   INTEGER, DIMENSION(:), ALLOCATABLE :: & 
     78      & n2dintsurf           !: Interpolation option for surface variables 
     79   REAL(wp), DIMENSION(:), ALLOCATABLE :: & 
     80      & zavglamscl, &        !: E/W diameter of averaging footprint for surface variables 
     81      & zavgphiscl           !: N/S diameter of averaging footprint for surface variables 
     82   LOGICAL, DIMENSION(:), ALLOCATABLE :: & 
     83      & lfpindegs, &         !: T=> surface obs footprint size specified in degrees, F=> in metres 
     84      & llnightav            !: Logical for calculating night-time averages 
     85 
     86   TYPE(obs_surf), PUBLIC, POINTER, DIMENSION(:) :: & 
     87      & surfdata, &          !: Initial surface data 
     88      & surfdataqc           !: Surface data after quality control 
     89   TYPE(obs_prof), PUBLIC, POINTER, DIMENSION(:) :: & 
     90      & profdata, &          !: Initial profile data 
     91      & profdataqc           !: Profile data after quality control 
    5992 
    6093   CHARACTER(len=6), PUBLIC, DIMENSION(:), ALLOCATABLE ::   cobstypesprof, cobstypessurf   !: Profile & surface obs types 
     
    95128         & cn_profbfiles, &      ! T/S profile input filenames 
    96129         & cn_sstfbfiles, &      ! Sea surface temperature input filenames 
     130         & cn_sssfbfiles, &      ! Sea surface salinity input filenames 
    97131         & cn_slafbfiles, &      ! Sea level anomaly input filenames 
    98132         & cn_sicfbfiles, &      ! Seaice concentration input filenames 
    99133         & cn_velfbfiles, &      ! Velocity profile input filenames 
    100          & cn_sstbias_files      ! SST bias input filenames 
     134         & cn_sstbiasfiles      ! SST bias input filenames 
    101135      CHARACTER(LEN=128) :: & 
    102136         & cn_altbiasfile        ! Altimeter bias input filename 
     
    109143      LOGICAL :: ln_sla          ! Logical switch for sea level anomalies  
    110144      LOGICAL :: ln_sst          ! Logical switch for sea surface temperature 
     145      LOGICAL :: ln_sss          ! Logical switch for sea surface salinity 
    111146      LOGICAL :: ln_sic          ! Logical switch for sea ice concentration 
    112147      LOGICAL :: ln_vel3d        ! Logical switch for velocity (u,v) obs 
     
    116151      LOGICAL :: ln_ignmis       ! Logical switch for ignoring missing files 
    117152      LOGICAL :: ln_s_at_t       ! Logical switch to compute model S at T obs 
     153      LOGICAL :: ln_bound_reject ! Logical to remove obs near boundaries in LAMs. 
    118154      LOGICAL :: llvar1          ! Logical for profile variable 1 
    119155      LOGICAL :: llvar2          ! Logical for profile variable 1 
    120       LOGICAL :: llnightav       ! Logical for calculating night-time averages 
    121156      LOGICAL, DIMENSION(jpmaxnfiles) :: lmask ! Used for finding number of sstbias files 
    122157 
     
    134169 
    135170      NAMELIST/namobs/ln_diaobs, ln_t3d, ln_s3d, ln_sla,              & 
    136          &            ln_sst, ln_sic, ln_vel3d,                       & 
    137          &            ln_altbias, ln_nea, ln_grid_global,             & 
    138          &            ln_grid_search_lookup,                          & 
    139          &            ln_ignmis, ln_s_at_t, ln_sstnight,              & 
     171         &            ln_sst, ln_sic, ln_sss, ln_vel3d,               & 
     172         &            ln_altbias, ln_sstbias, ln_nea,                 & 
     173         &            ln_grid_global, ln_grid_search_lookup,          & 
     174         &            ln_ignmis, ln_s_at_t, ln_bound_reject,          & 
     175         &            ln_sstnight,                                    & 
     176         &            ln_sla_fp_indegs, ln_sst_fp_indegs,             & 
     177         &            ln_sss_fp_indegs, ln_sic_fp_indegs,             & 
    140178         &            cn_profbfiles, cn_slafbfiles,                   & 
    141179         &            cn_sstfbfiles, cn_sicfbfiles,                   & 
    142          &            cn_velfbfiles, cn_altbiasfile,                  & 
     180         &            cn_velfbfiles, cn_sssfbfiles,                   & 
     181         &            cn_sstbiasfiles, cn_altbiasfile,                & 
    143182         &            cn_gridsearchfile, rn_gridsearchres,            & 
    144          &            rn_dobsini, rn_dobsend, nn_1dint, nn_2dint,     & 
     183         &            rn_dobsini, rn_dobsend,                         & 
     184         &            rn_sla_avglamscl, rn_sla_avgphiscl,             & 
     185         &            rn_sst_avglamscl, rn_sst_avgphiscl,             & 
     186         &            rn_sss_avglamscl, rn_sss_avgphiscl,             & 
     187         &            rn_sic_avglamscl, rn_sic_avgphiscl,             & 
     188         &            nn_1dint, nn_2dint,                             & 
     189         &            nn_2dint_sla, nn_2dint_sst,                     & 
     190         &            nn_2dint_sss, nn_2dint_sic,                     & 
    145191         &            nn_msshc, rn_mdtcorr, rn_mdtcutoff,             & 
    146          &            nn_profdavtypes, ln_sstbias, cn_sstbias_files 
     192         &            nn_profdavtypes 
    147193 
    148194      INTEGER :: jnumsstbias 
     
    157203      ! Read namelist parameters 
    158204      !----------------------------------------------------------------------- 
    159        
    160       !Initalise all values in namelist arrays 
    161       ALLOCATE(sstbias_type(jpmaxnfiles)) 
    162205      ! Some namelist arrays need initialising 
    163206      cn_profbfiles(:) = '' 
     
    166209      cn_sicfbfiles(:) = '' 
    167210      cn_velfbfiles(:) = '' 
    168       cn_sstbias_files(:) = '' 
     211      cn_sssfbfiles(:)    = '' 
     212      cn_sstbiasfiles(:) = '' 
    169213      nn_profdavtypes(:) = -1 
    170214 
     
    187231         RETURN 
    188232      ENDIF 
    189        
    190       !----------------------------------------------------------------------- 
    191       ! Set up list of observation types to be used 
    192       ! and the files associated with each type 
    193       !----------------------------------------------------------------------- 
    194  
    195       nproftypes = COUNT( (/ln_t3d .OR. ln_s3d, ln_vel3d /) ) 
    196       nsurftypes = COUNT( (/ln_sla, ln_sst, ln_sic /) ) 
    197  
    198       IF (ln_sstbias) THEN  
    199          lmask(:) = .FALSE.  
    200          WHERE (cn_sstbias_files(:) /= '') lmask(:) = .TRUE.  
    201          jnumsstbias = COUNT(lmask)  
    202          lmask(:) = .FALSE.  
    203       ENDIF       
    204  
    205       IF ( nproftypes == 0 .AND. nsurftypes == 0 ) THEN 
    206          IF(lwp) WRITE(numout,cform_war) 
    207          IF(lwp) WRITE(numout,*) ' ln_diaobs is set to true, but all obs operator logical flags', & 
    208             &                    ' ln_t3d, ln_s3d, ln_sla, ln_sst, ln_sic, ln_vel3d', & 
    209             &                    ' are set to .FALSE. so turning off calls to dia_obs' 
    210          nwarn = nwarn + 1 
    211          ln_diaobs = .FALSE. 
    212          RETURN 
    213       ENDIF 
    214  
    215       IF ( nproftypes > 0 ) THEN 
    216  
    217          ALLOCATE( cobstypesprof(nproftypes) ) 
    218          ALLOCATE( ifilesprof(nproftypes) ) 
    219          ALLOCATE( clproffiles(nproftypes,jpmaxnfiles) ) 
    220  
    221          jtype = 0 
    222          IF (ln_t3d .OR. ln_s3d) THEN 
    223             jtype = jtype + 1 
    224             clproffiles(jtype,:) = cn_profbfiles(:) 
    225             cobstypesprof(jtype) = 'prof  ' 
    226             ifilesprof(jtype) = 0 
    227             DO jfile = 1, jpmaxnfiles 
    228                IF ( trim(clproffiles(jtype,jfile)) /= '' ) & 
    229                   ifilesprof(jtype) = ifilesprof(jtype) + 1 
    230             END DO 
    231          ENDIF 
    232          IF (ln_vel3d) THEN 
    233             jtype = jtype + 1 
    234             clproffiles(jtype,:) = cn_velfbfiles(:) 
    235             cobstypesprof(jtype) = 'vel   ' 
    236             ifilesprof(jtype) = 0 
    237             DO jfile = 1, jpmaxnfiles 
    238                IF ( trim(clproffiles(jtype,jfile)) /= '' ) & 
    239                   ifilesprof(jtype) = ifilesprof(jtype) + 1 
    240             END DO 
    241          ENDIF 
    242  
    243       ENDIF 
    244  
    245       IF ( nsurftypes > 0 ) THEN 
    246  
    247          ALLOCATE( cobstypessurf(nsurftypes) ) 
    248          ALLOCATE( ifilessurf(nsurftypes) ) 
    249          ALLOCATE( clsurffiles(nsurftypes, jpmaxnfiles) ) 
    250  
    251          jtype = 0 
    252          IF (ln_sla) THEN 
    253             jtype = jtype + 1 
    254             clsurffiles(jtype,:) = cn_slafbfiles(:) 
    255             cobstypessurf(jtype) = 'sla   ' 
    256             ifilessurf(jtype) = 0 
    257             DO jfile = 1, jpmaxnfiles 
    258                IF ( trim(clsurffiles(jtype,jfile)) /= '' ) & 
    259                   ifilessurf(jtype) = ifilessurf(jtype) + 1 
    260             END DO 
    261          ENDIF 
    262          IF (ln_sst) THEN 
    263             jtype = jtype + 1 
    264             clsurffiles(jtype,:) = cn_sstfbfiles(:) 
    265             cobstypessurf(jtype) = 'sst   ' 
    266             ifilessurf(jtype) = 0 
    267             DO jfile = 1, jpmaxnfiles 
    268                IF ( trim(clsurffiles(jtype,jfile)) /= '' ) & 
    269                   ifilessurf(jtype) = ifilessurf(jtype) + 1 
    270             END DO 
    271          ENDIF 
    272 #if defined key_lim3 
    273          IF (ln_sic) THEN 
    274             jtype = jtype + 1 
    275             clsurffiles(jtype,:) = cn_sicfbfiles(:) 
    276             cobstypessurf(jtype) = 'sic   ' 
    277             ifilessurf(jtype) = 0 
    278             DO jfile = 1, jpmaxnfiles 
    279                IF ( trim(clsurffiles(jtype,jfile)) /= '' ) & 
    280                   ifilessurf(jtype) = ifilessurf(jtype) + 1 
    281             END DO 
    282          ENDIF 
    283 #endif 
    284  
    285       ENDIF 
    286  
    287       !Write namelist settings to stdout 
     233 
    288234      IF(lwp) THEN 
    289235         WRITE(numout,*) 
     
    297243         WRITE(numout,*) '             Logical switch for Sea Ice observations                  ln_sic = ', ln_sic 
    298244         WRITE(numout,*) '             Logical switch for velocity observations               ln_vel3d = ', ln_vel3d 
    299          WRITE(numout,*) '             Global distribution of observations              ln_grid_global = ',ln_grid_global 
    300          WRITE(numout,*) '             Logical switch for SST bias correction         ln_sstbias = ', ln_sstbias  
    301          WRITE(numout,*) '             Logical switch for obs grid search lookup ln_grid_search_lookup = ',ln_grid_search_lookup 
     245         WRITE(numout,*) '             Logical switch for SSS observations                      ln_sss = ', ln_sss 
     246         WRITE(numout,*) '             Global distribution of observations              ln_grid_global = ', ln_grid_global 
     247         WRITE(numout,*) '             Logical switch for obs grid search lookup ln_grid_search_lookup = ', ln_grid_search_lookup 
    302248         IF (ln_grid_search_lookup) & 
    303249            WRITE(numout,*) '             Grid search lookup file header                cn_gridsearchfile = ', cn_gridsearchfile 
     
    307253         WRITE(numout,*) '             Type of horizontal interpolation method                nn_2dint = ', nn_2dint 
    308254         WRITE(numout,*) '             Rejection of observations near land switch               ln_nea = ', ln_nea 
     255         WRITE(numout,*) '             Rejection of obs near open bdys                 ln_bound_reject = ', ln_bound_reject 
    309256         WRITE(numout,*) '             MSSH correction scheme                                 nn_msshc = ', nn_msshc 
    310257         WRITE(numout,*) '             MDT  correction                                      rn_mdtcorr = ', rn_mdtcorr 
    311258         WRITE(numout,*) '             MDT cutoff for computed correction                 rn_mdtcutoff = ', rn_mdtcutoff 
    312259         WRITE(numout,*) '             Logical switch for alt bias                          ln_altbias = ', ln_altbias 
     260         WRITE(numout,*) '             Logical switch for sst bias                          ln_sstbias = ', ln_sstbias 
    313261         WRITE(numout,*) '             Logical switch for ignoring missing files             ln_ignmis = ', ln_ignmis 
    314262         WRITE(numout,*) '             Daily average types                             nn_profdavtypes = ', nn_profdavtypes 
    315263         WRITE(numout,*) '             Logical switch for night-time SST obs               ln_sstnight = ', ln_sstnight 
    316          WRITE(numout,*) '          Number of profile obs types: ',nproftypes 
    317  
    318          IF ( nproftypes > 0 ) THEN 
    319             DO jtype = 1, nproftypes 
    320                DO jfile = 1, ifilesprof(jtype) 
    321                   WRITE(numout,'(1X,2A)') '             '//cobstypesprof(jtype)//' input observation file names  = ', & 
    322                      TRIM(clproffiles(jtype,jfile)) 
    323                END DO 
    324             END DO 
     264      ENDIF 
     265      !----------------------------------------------------------------------- 
     266      ! Set up list of observation types to be used 
     267      ! and the files associated with each type 
     268      !----------------------------------------------------------------------- 
     269 
     270      nproftypes = COUNT( (/ln_t3d .OR. ln_s3d, ln_vel3d /) ) 
     271      nsurftypes = COUNT( (/ln_sla, ln_sst, ln_sic, ln_sss /) ) 
     272 
     273      IF (ln_sstbias) THEN  
     274         lmask(:) = .FALSE.  
     275         WHERE (cn_sstbiasfiles(:) /= '') lmask(:) = .TRUE.  
     276         jnumsstbias = COUNT(lmask)  
     277         lmask(:) = .FALSE.  
     278      ENDIF       
     279 
     280      IF ( nproftypes == 0 .AND. nsurftypes == 0 ) THEN 
     281         IF(lwp) WRITE(numout,cform_war) 
     282         IF(lwp) WRITE(numout,*) ' ln_diaobs is set to true, but all obs operator logical flags', & 
     283            &                    ' ln_t3d, ln_s3d, ln_sla, ln_sst, ln_sic, ln_vel3d', & 
     284            &                    ' are set to .FALSE. so turning off calls to dia_obs' 
     285         nwarn = nwarn + 1 
     286         ln_diaobs = .FALSE. 
     287         RETURN 
     288      ENDIF 
     289 
     290      IF ( nproftypes > 0 ) THEN 
     291 
     292         ALLOCATE( cobstypesprof(nproftypes) ) 
     293         ALLOCATE( ifilesprof(nproftypes) ) 
     294         ALLOCATE( clproffiles(nproftypes,jpmaxnfiles) ) 
     295 
     296         jtype = 0 
     297         IF (ln_t3d .OR. ln_s3d) THEN 
     298            jtype = jtype + 1 
     299            CALL obs_settypefiles( nproftypes, jpmaxnfiles, jtype, 'prof  ', & 
     300               &                   cn_profbfiles, ifilesprof, cobstypesprof, clproffiles ) 
    325301         ENDIF 
    326  
    327          WRITE(numout,*)'          Number of surface obs types: ',nsurftypes 
    328          IF ( nsurftypes > 0 ) THEN 
    329             DO jtype = 1, nsurftypes 
    330                DO jfile = 1, ifilessurf(jtype) 
    331                   WRITE(numout,'(1X,2A)') '             '//cobstypessurf(jtype)//' input observation file names  = ', & 
    332                      TRIM(clsurffiles(jtype,jfile)) 
    333                END DO 
    334             END DO 
     302         IF (ln_vel3d) THEN 
     303            jtype = jtype + 1 
     304            CALL obs_settypefiles( nproftypes, jpmaxnfiles, jtype, 'vel   ', & 
     305               &                   cn_velfbfiles, ifilesprof, cobstypesprof, clproffiles ) 
    335306         ENDIF 
    336          WRITE(numout,*) '~~~~~~~~~~~~' 
    337  
    338       ENDIF 
     307 
     308      ENDIF 
     309 
     310      IF ( nsurftypes > 0 ) THEN 
     311 
     312         ALLOCATE( cobstypessurf(nsurftypes) ) 
     313         ALLOCATE( ifilessurf(nsurftypes) ) 
     314         ALLOCATE( clsurffiles(nsurftypes, jpmaxnfiles) ) 
     315         ALLOCATE(n2dintsurf(nsurftypes)) 
     316         ALLOCATE(zavglamscl(nsurftypes)) 
     317         ALLOCATE(zavgphiscl(nsurftypes)) 
     318         ALLOCATE(lfpindegs(nsurftypes)) 
     319         ALLOCATE(llnightav(nsurftypes)) 
     320 
     321         jtype = 0 
     322         IF (ln_sla) THEN 
     323            jtype = jtype + 1 
     324            CALL obs_settypefiles( nsurftypes, jpmaxnfiles, jtype, 'sla   ', & 
     325               &                   cn_slafbfiles, ifilessurf, cobstypessurf, clsurffiles ) 
     326            CALL obs_setinterpopts( nsurftypes, jtype, 'sla   ',      & 
     327               &                  nn_2dint, nn_2dint_sla,             & 
     328               &                  rn_sla_avglamscl, rn_sla_avgphiscl, & 
     329               &                  ln_sla_fp_indegs, .FALSE.,          & 
     330               &                  n2dintsurf, zavglamscl, zavgphiscl, & 
     331               &                  lfpindegs, llnightav ) 
     332         ENDIF 
     333         IF (ln_sst) THEN 
     334            jtype = jtype + 1 
     335            CALL obs_settypefiles( nsurftypes, jpmaxnfiles, jtype, 'sst   ', & 
     336               &                   cn_sstfbfiles, ifilessurf, cobstypessurf, clsurffiles ) 
     337            CALL obs_setinterpopts( nsurftypes, jtype, 'sst   ',      & 
     338               &                  nn_2dint, nn_2dint_sst,             & 
     339               &                  rn_sst_avglamscl, rn_sst_avgphiscl, & 
     340               &                  ln_sst_fp_indegs, ln_sstnight,      & 
     341               &                  n2dintsurf, zavglamscl, zavgphiscl, & 
     342               &                  lfpindegs, llnightav ) 
     343         ENDIF 
     344#if defined key_lim3 || defined key_cice 
     345         IF (ln_sic) THEN 
     346            jtype = jtype + 1 
     347            CALL obs_settypefiles( nsurftypes, jpmaxnfiles, jtype, 'sic   ', & 
     348               &                   cn_sicfbfiles, ifilessurf, cobstypessurf, clsurffiles ) 
     349            CALL obs_setinterpopts( nsurftypes, jtype, 'sic   ',      & 
     350               &                  nn_2dint, nn_2dint_sic,             & 
     351               &                  rn_sic_avglamscl, rn_sic_avgphiscl, & 
     352               &                  ln_sic_fp_indegs, .FALSE.,          & 
     353               &                  n2dintsurf, zavglamscl, zavgphiscl, & 
     354               &                  lfpindegs, llnightav ) 
     355         ENDIF 
     356#endif 
     357         IF (ln_sss) THEN 
     358            jtype = jtype + 1 
     359            CALL obs_settypefiles( nsurftypes, jpmaxnfiles, jtype, 'sss   ', & 
     360               &                   cn_sssfbfiles, ifilessurf, cobstypessurf, clsurffiles ) 
     361            CALL obs_setinterpopts( nsurftypes, jtype, 'sss   ',      & 
     362               &                  nn_2dint, nn_2dint_sss,             & 
     363               &                  rn_sss_avglamscl, rn_sss_avgphiscl, & 
     364               &                  ln_sss_fp_indegs, .FALSE.,          & 
     365               &                  n2dintsurf, zavglamscl, zavgphiscl, & 
     366               &                  lfpindegs, llnightav ) 
     367         ENDIF 
     368 
     369      ENDIF 
     370 
     371 
    339372 
    340373      !----------------------------------------------------------------------- 
     
    356389      ENDIF 
    357390 
    358       IF ( ( nn_2dint < 0 ) .OR. ( nn_2dint > 4 ) ) THEN 
     391      IF ( ( nn_2dint < 0 ) .OR. ( nn_2dint > 6 ) ) THEN 
    359392         CALL ctl_stop(' Choice of horizontal (2D) interpolation method', & 
    360393            &                    ' is not available') 
     
    421454               &               jpi, jpj, jpk, & 
    422455               &               zmask1, zglam1, zgphi1, zmask2, zglam2, zgphi2,  & 
    423                &               ln_nea, kdailyavtypes = nn_profdavtypes ) 
     456               &               ln_nea, ln_bound_reject, & 
     457               &               kdailyavtypes = nn_profdavtypes ) 
    424458 
    425459         END DO 
     
    440474            nvarssurf(jtype) = 1 
    441475            nextrsurf(jtype) = 0 
    442             llnightav = .FALSE. 
     476            llnightav(jtype) = .FALSE. 
    443477            IF ( TRIM(cobstypessurf(jtype)) == 'sla' ) nextrsurf(jtype) = 2 
    444             IF ( TRIM(cobstypessurf(jtype)) == 'sst' ) llnightav = ln_sstnight 
     478            IF ( TRIM(cobstypessurf(jtype)) == 'sst' ) llnightav(jtype) = ln_sstnight 
    445479 
    446480            !Read in surface obs types 
     
    448482               &               clsurffiles(jtype,1:ifilessurf(jtype)), & 
    449483               &               nvarssurf(jtype), nextrsurf(jtype), nitend-nit000+2, & 
    450                &               rn_dobsini, rn_dobsend, ln_ignmis, .FALSE., llnightav ) 
    451           
    452           
    453             CALL obs_pre_surf( surfdata(jtype), surfdataqc(jtype), ln_nea ) 
     484               &               rn_dobsini, rn_dobsend, ln_ignmis, .FALSE., llnightav(jtype) ) 
     485 
     486            CALL obs_pre_surf( surfdata(jtype), surfdataqc(jtype), ln_nea, ln_bound_reject ) 
    454487 
    455488            IF ( TRIM(cobstypessurf(jtype)) == 'sla' ) THEN 
    456                CALL obs_rea_mdt( surfdataqc(jtype), nn_2dint ) 
    457                IF ( ln_altbias ) CALL obs_rea_altbias ( surfdataqc(jtype), nn_2dint, cn_altbiasfile ) 
     489               CALL obs_rea_mdt( surfdataqc(jtype), n2dintsurf(jtype) ) 
     490               IF ( ln_altbias ) & 
     491                  & CALL obs_rea_altbias ( surfdataqc(jtype), n2dintsurf(jtype), cn_altbiasfile ) 
    458492            ENDIF 
     493 
    459494            IF ( TRIM(cobstypessurf(jtype)) == 'sst' .AND. ln_sstbias ) THEN 
    460                !Read in bias field and correct SST. 
    461                IF ( jnumsstbias == 0 ) CALL ctl_stop("ln_sstbias set,"// & 
    462                                                      "  but no bias"// & 
    463                                                      " files to read in")    
    464                   CALL obs_app_sstbias( surfdataqc(jtype), nn_2dint, & 
    465                                         jnumsstbias, cn_sstbias_files(1:jnumsstbias) ) 
     495               jnumsstbias = 0 
     496               DO jfile = 1, jpmaxnfiles 
     497                  IF ( TRIM(cn_sstbiasfiles(jfile)) /= '' ) & 
     498                     &  jnumsstbias = jnumsstbias + 1 
     499               END DO 
     500               IF ( jnumsstbias == 0 ) THEN 
     501                  CALL ctl_stop("ln_sstbias set but no bias files to read in")     
     502               ENDIF 
     503 
     504               CALL obs_app_sstbias( surfdataqc(jtype), n2dintsurf(jtype), &  
     505                  &                  jnumsstbias, cn_sstbiasfiles(1:jnumsstbias) )  
     506 
    466507            ENDIF 
    467508         END DO 
     
    512553      USE ice    , ONLY : at_i                ! LIM3 Ice model variables 
    513554#endif 
     555#if defined key_cice 
     556      USE sbc_oce, ONLY : fr_i     ! ice fraction 
     557#endif 
     558 
    514559      IMPLICIT NONE 
    515560 
     
    528573         & zprofmask2              ! Mask associated with zprofvar2 
    529574      REAL(wp), POINTER, DIMENSION(:,:) :: & 
    530          & zsurfvar                ! Model values equivalent to surface ob. 
     575         & zsurfvar, &             ! Model values equivalent to surface ob. 
     576         & zsurfmask               ! Mask associated with surface variable 
    531577      REAL(wp), POINTER, DIMENSION(:,:) :: & 
    532578         & zglam1,    &            ! Model longitudes for prof variable 1 
     
    534580         & zgphi1,    &            ! Model latitudes for prof variable 1 
    535581         & zgphi2                  ! Model latitudes for prof variable 2 
    536 #if ! defined key_lim3 
    537       REAL(wp), POINTER, DIMENSION(:,:) :: at_i 
    538 #endif 
    539       LOGICAL :: llnightav        ! Logical for calculating night-time average 
    540582 
    541583      !Allocate local work arrays 
     
    545587      CALL wrk_alloc( jpi, jpj, jpk, zprofmask2 ) 
    546588      CALL wrk_alloc( jpi, jpj, zsurfvar ) 
     589      CALL wrk_alloc( jpi, jpj, zsurfmask ) 
    547590      CALL wrk_alloc( jpi, jpj, zglam1 ) 
    548591      CALL wrk_alloc( jpi, jpj, zglam2 ) 
    549592      CALL wrk_alloc( jpi, jpj, zgphi1 ) 
    550593      CALL wrk_alloc( jpi, jpj, zgphi2 ) 
    551 #if ! defined key_lim3 
    552       CALL wrk_alloc(jpi,jpj,at_i)  
    553 #endif 
    554594      !----------------------------------------------------------------------- 
    555595 
     
    562602      idaystp = NINT( rday / rdt ) 
    563603 
    564       !----------------------------------------------------------------------- 
    565       ! No LIM => at_i == 0.0_wp 
    566       !----------------------------------------------------------------------- 
    567 #if ! defined key_lim3 
    568       at_i(:,:) = 0.0_wp 
    569 #endif 
    570604      !----------------------------------------------------------------------- 
    571605      ! Call the profile and surface observation operators 
     
    595629               zgphi1(:,:) = gphiu(:,:) 
    596630               zgphi2(:,:) = gphiv(:,:) 
     631            CASE DEFAULT 
     632               CALL ctl_stop( 'Unknown profile observation type '//TRIM(cobstypesprof(jtype))//' in dia_obs' ) 
    597633            END SELECT 
    598634 
    599             IF( ln_zco .OR. ln_zps ) THEN  
    600                CALL obs_prof_opt( profdataqc(jtype), kstp, jpi, jpj, jpk,  & 
    601                   &               nit000, idaystp,                         & 
    602                   &               zprofvar1, zprofvar2,                    & 
    603                   &               gdept_1d, zprofmask1, zprofmask2,        & 
    604                   &               zglam1, zglam2, zgphi1, zgphi2,          & 
    605                   &               nn_1dint, nn_2dint,                      & 
    606                   &               kdailyavtypes = nn_profdavtypes ) 
    607             ELSE IF(TRIM(cobstypesprof(jtype)) == 'prof') THEN 
    608                !TG - THIS NEEDS MODIFICATION TO MATCH SIMPLIFICATION 
    609                CALL obs_pro_sco_opt( profdataqc(jtype),                    &  
    610                   &              kstp, jpi, jpj, jpk, nit000, idaystp,   &  
    611                   &              zprofvar1, zprofvar2,                   &  
    612                   &              gdept_n(:,:,:), gdepw_n(:,:,:),           & 
    613                   &              tmask, nn_1dint, nn_2dint,              &  
    614                   &              kdailyavtypes = nn_profdavtypes )  
    615             ELSE 
    616                CALL ctl_stop('DIA_OBS: Generalised vertical interpolation not'// & 
    617                          'yet working for velocity data (turn off velocity observations') 
    618             ENDIF 
     635            CALL obs_prof_opt( profdataqc(jtype), kstp, jpi, jpj, jpk,  & 
     636               &               nit000, idaystp,                         & 
     637               &               zprofvar1, zprofvar2,                    & 
     638               &               gdept_n(:,:,:), gdepw_n(:,:,:),            &  
     639               &               zprofmask1, zprofmask2,                  & 
     640               &               zglam1, zglam2, zgphi1, zgphi2,          & 
     641               &               nn_1dint, nn_2dint,                      & 
     642               &               kdailyavtypes = nn_profdavtypes ) 
    619643 
    620644         END DO 
     
    625649 
    626650         DO jtype = 1, nsurftypes 
     651 
     652            !Defaults which might be changed 
     653            zsurfmask(:,:) = tmask(:,:,1) 
    627654 
    628655            SELECT CASE ( TRIM(cobstypessurf(jtype)) ) 
    629656            CASE('sst') 
    630657               zsurfvar(:,:) = tsn(:,:,1,jp_tem) 
    631                llnightav = ln_sstnight 
    632658            CASE('sla') 
    633659               zsurfvar(:,:) = sshn(:,:) 
    634                llnightav = .FALSE. 
    635 #if defined key_lim3 
     660            CASE('sss') 
     661               zsurfvar(:,:) = tsn(:,:,1,jp_sal) 
    636662            CASE('sic') 
    637663               IF ( kstp == 0 ) THEN 
     
    646672                  CYCLE 
    647673               ELSE 
    648                   zsurfvar(:,:) = at_i(:,:) 
     674#if defined key_cice 
     675                  zsurfvar(:,:) = fr_i(:,:) 
     676#elif defined key_lim2 || defined key_lim3 
     677                  zsurfvar(:,:) = 1._wp - frld(:,:) 
     678#else 
     679               CALL ctl_stop( ' Trying to run sea-ice observation operator', & 
     680                  &           ' but no sea-ice model appears to have been defined' ) 
     681#endif 
    649682               ENDIF 
    650683 
    651                llnightav = .FALSE. 
    652 #endif 
    653684            END SELECT 
    654685 
    655686            CALL obs_surf_opt( surfdataqc(jtype), kstp, jpi, jpj,       & 
    656                &               nit000, idaystp, zsurfvar, tmask(:,:,1), & 
    657                &               nn_2dint, llnightav ) 
     687               &               nit000, idaystp, zsurfvar, zsurfmask,    & 
     688               &               n2dintsurf(jtype), llnightav(jtype),     & 
     689               &               zavglamscl(jtype), zavgphiscl(jtype),     & 
     690               &               lfpindegs(jtype) ) 
    658691 
    659692         END DO 
     
    666699      CALL wrk_dealloc( jpi, jpj, jpk, zprofmask2 ) 
    667700      CALL wrk_dealloc( jpi, jpj, zsurfvar ) 
     701      CALL wrk_dealloc( jpi, jpj, zsurfmask ) 
    668702      CALL wrk_dealloc( jpi, jpj, zglam1 ) 
    669703      CALL wrk_dealloc( jpi, jpj, zglam2 ) 
    670704      CALL wrk_dealloc( jpi, jpj, zgphi1 ) 
    671705      CALL wrk_dealloc( jpi, jpj, zgphi2 ) 
    672 #if ! defined key_lim3 
    673       CALL wrk_dealloc(jpi,jpj,at_i) 
    674 #endif 
    675706 
    676707   END SUBROUTINE dia_obs 
     
    789820 
    790821      IF ( nsurftypes > 0 ) & 
    791          &   DEALLOCATE( cobstypessurf, surfdata, surfdataqc, nvarssurf, nextrsurf ) 
     822         &   DEALLOCATE( cobstypessurf, surfdata, surfdataqc, nvarssurf, nextrsurf, & 
     823         &               n2dintsurf, zavglamscl, zavgphiscl, lfpindegs, llnightav ) 
    792824 
    793825   END SUBROUTINE dia_obs_dealloc 
     
    938970   END SUBROUTINE fin_date 
    939971    
     972    SUBROUTINE obs_settypefiles( ntypes, jpmaxnfiles, jtype, ctypein, & 
     973       &                         cfilestype, ifiles, cobstypes, cfiles ) 
     974 
     975    INTEGER, INTENT(IN) :: ntypes      ! Total number of obs types 
     976    INTEGER, INTENT(IN) :: jpmaxnfiles ! Maximum number of files allowed for each type 
     977    INTEGER, INTENT(IN) :: jtype       ! Index of the current type of obs 
     978    INTEGER, DIMENSION(ntypes), INTENT(INOUT) :: & 
     979       &                   ifiles      ! Out appended number of files for this type 
     980 
     981    CHARACTER(len=6), INTENT(IN) :: ctypein  
     982    CHARACTER(len=128), DIMENSION(jpmaxnfiles), INTENT(IN) :: & 
     983       &                   cfilestype  ! In list of files for this obs type 
     984    CHARACTER(len=6), DIMENSION(ntypes), INTENT(INOUT) :: & 
     985       &                   cobstypes   ! Out appended list of obs types 
     986    CHARACTER(len=128), DIMENSION(ntypes, jpmaxnfiles), INTENT(INOUT) :: & 
     987       &                   cfiles      ! Out appended list of files for all types 
     988 
     989    !Local variables 
     990    INTEGER :: jfile 
     991 
     992    cfiles(jtype,:) = cfilestype(:) 
     993    cobstypes(jtype) = ctypein 
     994    ifiles(jtype) = 0 
     995    DO jfile = 1, jpmaxnfiles 
     996       IF ( trim(cfiles(jtype,jfile)) /= '' ) & 
     997                 ifiles(jtype) = ifiles(jtype) + 1 
     998    END DO 
     999 
     1000    IF ( ifiles(jtype) == 0 ) THEN 
     1001         CALL ctl_stop( 'Logical for observation type '//TRIM(ctypein)//   & 
     1002            &           ' set to true but no files available to read' ) 
     1003    ENDIF 
     1004 
     1005    IF(lwp) THEN     
     1006       WRITE(numout,*) '             '//cobstypes(jtype)//' input observation file names:' 
     1007       DO jfile = 1, ifiles(jtype) 
     1008          WRITE(numout,*) '                '//TRIM(cfiles(jtype,jfile)) 
     1009       END DO 
     1010    ENDIF 
     1011 
     1012    END SUBROUTINE obs_settypefiles 
     1013 
     1014    SUBROUTINE obs_setinterpopts( ntypes, jtype, ctypein,             & 
     1015               &                  n2dint_default, n2dint_type,        & 
     1016               &                  zavglamscl_type, zavgphiscl_type,   & 
     1017               &                  lfp_indegs_type, lavnight_type,     & 
     1018               &                  n2dint, zavglamscl, zavgphiscl,     & 
     1019               &                  lfpindegs, lavnight ) 
     1020 
     1021    INTEGER, INTENT(IN)  :: ntypes             ! Total number of obs types 
     1022    INTEGER, INTENT(IN)  :: jtype              ! Index of the current type of obs 
     1023    INTEGER, INTENT(IN)  :: n2dint_default     ! Default option for interpolation type 
     1024    INTEGER, INTENT(IN)  :: n2dint_type        ! Option for interpolation type 
     1025    REAL(wp), INTENT(IN) :: & 
     1026       &                    zavglamscl_type, & !E/W diameter of obs footprint for this type 
     1027       &                    zavgphiscl_type    !N/S diameter of obs footprint for this type 
     1028    LOGICAL, INTENT(IN)  :: lfp_indegs_type    !T=> footprint in degrees, F=> in metres 
     1029    LOGICAL, INTENT(IN)  :: lavnight_type      !T=> obs represent night time average 
     1030    CHARACTER(len=6), INTENT(IN) :: ctypein  
     1031 
     1032    INTEGER, DIMENSION(ntypes), INTENT(INOUT) :: & 
     1033       &                    n2dint  
     1034    REAL(wp), DIMENSION(ntypes), INTENT(INOUT) :: & 
     1035       &                    zavglamscl, zavgphiscl 
     1036    LOGICAL, DIMENSION(ntypes), INTENT(INOUT) :: & 
     1037       &                    lfpindegs, lavnight 
     1038 
     1039    lavnight(jtype) = lavnight_type 
     1040 
     1041    IF ( (n2dint_type >= 1) .AND. (n2dint_type <= 6) ) THEN 
     1042       n2dint(jtype) = n2dint_type 
     1043    ELSE 
     1044       n2dint(jtype) = n2dint_default 
     1045    ENDIF 
     1046 
     1047    ! For averaging observation footprints set options for size of footprint  
     1048    IF ( (n2dint(jtype) > 4) .AND. (n2dint(jtype) <= 6) ) THEN 
     1049       IF ( zavglamscl_type > 0._wp ) THEN 
     1050          zavglamscl(jtype) = zavglamscl_type 
     1051       ELSE 
     1052          CALL ctl_stop( 'Incorrect value set for averaging footprint '// & 
     1053                         'scale (zavglamscl) for observation type '//TRIM(ctypein) )       
     1054       ENDIF 
     1055 
     1056       IF ( zavgphiscl_type > 0._wp ) THEN 
     1057          zavgphiscl(jtype) = zavgphiscl_type 
     1058       ELSE 
     1059          CALL ctl_stop( 'Incorrect value set for averaging footprint '// & 
     1060                         'scale (zavgphiscl) for observation type '//TRIM(ctypein) )       
     1061       ENDIF 
     1062 
     1063       lfpindegs(jtype) = lfp_indegs_type  
     1064 
     1065    ENDIF 
     1066 
     1067    ! Write out info  
     1068    IF(lwp) THEN 
     1069       IF ( n2dint(jtype) <= 4 ) THEN 
     1070          WRITE(numout,*) '             '//TRIM(ctypein)// & 
     1071             &            ' model counterparts will be interpolated horizontally' 
     1072       ELSE IF ( n2dint(jtype) <= 6 ) THEN 
     1073          WRITE(numout,*) '             '//TRIM(ctypein)// & 
     1074             &            ' model counterparts will be averaged horizontally' 
     1075          WRITE(numout,*) '             '//'    with E/W scale: ',zavglamscl(jtype) 
     1076          WRITE(numout,*) '             '//'    with N/S scale: ',zavgphiscl(jtype) 
     1077          IF ( lfpindegs(jtype) ) THEN 
     1078              WRITE(numout,*) '             '//'    (in degrees)' 
     1079          ELSE 
     1080              WRITE(numout,*) '             '//'    (in metres)' 
     1081          ENDIF 
     1082       ENDIF 
     1083    ENDIF 
     1084 
     1085    END SUBROUTINE obs_setinterpopts 
     1086 
    9401087END MODULE diaobs 
Note: See TracChangeset for help on using the changeset viewer.