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 6854 – NEMO

Changeset 6854


Ignore:
Timestamp:
2016-08-08T12:26:45+02:00 (8 years ago)
Author:
dford
Message:

Initial implementation of observation operator for LogChl?.

Location:
branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS
Files:
4 edited
4 copied

Legend:

Unmodified
Added
Removed
  • branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90

    r6406 r6854  
    2121   USE par_oce 
    2222   USE dom_oce                  ! Ocean space and time domain variables 
     23   USE obs_const, ONLY: obfillflt ! Fill value 
    2324   USE obs_fbm, ONLY: ln_cl4    ! Class 4 diagnostic switch 
    2425   USE obs_read_prof            ! Reading and allocation of observations (Coriolis) 
     
    2930   USE obs_read_seaice          ! Reading and allocation of Sea Ice observations   
    3031   USE obs_read_vel             ! Reading and allocation of velocity component observations 
     32   USE obs_read_logchl          ! Reading and allocation of logchl observations 
    3133   USE obs_prep                 ! Preparation of obs. (grid search etc). 
    3234   USE obs_oper                 ! Observation operators 
     
    4042   USE obs_sst                  ! SST data storage 
    4143   USE obs_seaice               ! Sea Ice data storage 
     44   USE obs_logchl               ! logchl data storage 
    4245   USE obs_types                ! Definitions for observation types 
    4346   USE mpp_map                  ! MPP mapping 
     
    8184   LOGICAL, PUBLIC :: ln_velhradcp   !: Logical switch for raw high freq netCDF ADCP vel. data  
    8285   LOGICAL, PUBLIC :: ln_velfb       !: Logical switch for velocities from feedback files 
     86   LOGICAL, PUBLIC :: ln_logchl      !: Logical switch for log10(chlorophyll) 
     87   LOGICAL, PUBLIC :: ln_logchlfb    !: Logical switch for logchl from feedback files 
    8388   LOGICAL, PUBLIC :: ln_ssh         !: Logical switch for sea surface height 
    8489   LOGICAL, PUBLIC :: ln_sss         !: Logical switch for sea surface salinity 
     
    164169      CHARACTER(len=128) :: velhradcpfiles(MaxNumFiles) 
    165170      CHARACTER(len=128) :: velfbfiles(MaxNumFiles) 
     171      CHARACTER(len=128) :: logchlfiles(MaxNumFiles) 
     172      CHARACTER(len=128) :: logchlfbfiles(MaxNumFiles) 
    166173      CHARACTER(LEN=128) :: reysstname 
    167174      CHARACTER(LEN=12)  :: reysstfmt 
     
    188195         &            ln_velhradcp, velhradcpfiles,                   & 
    189196         &            ln_velfb, velfbfiles, ln_velfb_av,              & 
     197         &            ln_logchl, ln_logchlfb,                         & 
     198         &            logchlfiles, logchlfbfiles,                     & 
    190199         &            ln_profb_enatim, ln_ignmis, ln_cl4,             & 
    191200         &            ln_sstbias, sstbias_files 
     
    209218      INTEGER :: jnumvelhradcp    
    210219      INTEGER :: jnumvelfb 
     220      INTEGER :: jnumlogchl 
     221      INTEGER :: jnumlogchlfb 
    211222      INTEGER :: ji 
    212223      INTEGER :: jset 
     
    217228      ! Read namelist parameters 
    218229      !----------------------------------------------------------------------- 
     230 
     231      ln_logchl   = .FALSE. 
     232      ln_logchlfb = .FALSE. 
    219233       
    220234      !Initalise all values in namelist arrays 
     
    237251      velcurfiles(:) = '' 
    238252      veladcpfiles(:) = '' 
     253      logchlfiles(:) = '' 
     254      logchlfbfiles(:) = '' 
    239255      sstbias_files(:) = '' 
    240256      endailyavtypes(:) = -1 
     
    335351         jnumvelfb = COUNT(lmask) 
    336352         lmask(:) = .FALSE. 
     353      ENDIF 
     354      IF (ln_logchl) THEN 
     355         lmask(:) = .FALSE. 
     356         WHERE (logchlfiles(:) /= '') lmask(:) = .TRUE. 
     357         jnumlogchl = COUNT(lmask) 
     358      ENDIF 
     359      IF (ln_logchlfb) THEN 
     360         lmask(:) = .FALSE. 
     361         WHERE (logchlfbfiles(:) /= '') lmask(:) = .TRUE. 
     362         jnumlogchlfb = COUNT(lmask) 
    337363      ENDIF 
    338364       
     
    366392         WRITE(numout,*) '             Logical switch for velocity high freq. ADCP  ln_velhradcp = ', ln_velhradcp 
    367393         WRITE(numout,*) '             Logical switch for feedback velocity data        ln_velfb = ', ln_velfb 
     394         WRITE(numout,*) '             Logical switch for logchl observations          ln_logchl = ', ln_logchl 
     395         WRITE(numout,*) '             Logical switch for feedback logchl data       ln_logchlfb = ', ln_logchlfb 
    368396         WRITE(numout,*) '             Global distribtion of observations         ln_grid_global = ',ln_grid_global 
    369397         WRITE(numout,*) & 
     
    462490                     TRIM(velfbfiles(ji)) 
    463491               ENDIF 
     492            END DO 
     493         ENDIF 
     494         IF (ln_logchl) THEN 
     495            DO ji = 1, jnumlogchl 
     496               WRITE(numout,'(1X,2A)') '             logchl input observation file name        logchlfiles = ', & 
     497                  TRIM(logchlfiles(ji)) 
     498            END DO 
     499         ENDIF 
     500         IF (ln_logchlfb) THEN 
     501            DO ji = 1, jnumlogchlfb 
     502               WRITE(numout,'(1X,2A)') '        Feedback logchl input observation file name  logchlfbfiles = ', & 
     503                  TRIM(logchlfbfiles(ji)) 
    464504            END DO 
    465505         ENDIF 
     
    498538         & ( .NOT. ln_vel3d ).AND.                                         & 
    499539         & ( .NOT. ln_ssh ).AND.( .NOT. ln_sst ).AND.( .NOT. ln_sss ).AND. & 
    500          & ( .NOT. ln_seaice ).AND.( .NOT. ln_vel3d ) ) THEN 
     540         & ( .NOT. ln_seaice ).AND.( .NOT. ln_vel3d ).AND.( .NOT. ln_logchl ) ) THEN 
    501541         IF(lwp) WRITE(numout,cform_war) 
    502542         IF(lwp) WRITE(numout,*) ' key_diaobs is activated but logical flags', & 
    503             &                    ' ln_t3d, ln_s3d, ln_sla, ln_ssh, ln_sst, ln_sss, ln_seaice, ln_vel3d are all set to .FALSE.' 
     543            &                    ' ln_t3d, ln_s3d, ln_sla, ln_ssh, ln_sst, ln_sss, ln_seaice, ln_vel3d,', & 
     544            &                    ' ln_logchl are all set to .FALSE.' 
    504545         nwarn = nwarn + 1 
    505546      ENDIF 
     
    9991040 
    10001041      ENDIF 
     1042 
     1043      !  - log10(chlorophyll) 
     1044       
     1045      IF ( ln_logchl ) THEN 
     1046 
     1047         ! Set the number of variables for logchl to 1 
     1048         nlogchlvars = 1 
     1049 
     1050         ! Set the number of extra variables for logchl to 0 
     1051         nlogchlextr = 0 
     1052          
     1053         IF ( ln_logchlfb ) THEN 
     1054            nlogchlsets = jnumlogchlfb 
     1055         ELSE 
     1056            nlogchlsets = 1 
     1057         ENDIF 
     1058 
     1059         ALLOCATE(logchldata(nlogchlsets)) 
     1060         ALLOCATE(logchldatqc(nlogchlsets)) 
     1061         logchldata(:)%nsurf=0 
     1062         logchldatqc(:)%nsurf=0 
     1063 
     1064         nlogchlsets = 0 
     1065 
     1066         IF ( ln_logchlfb ) THEN             ! Feedback file format 
     1067 
     1068            DO jset = 1, jnumlogchlfb 
     1069             
     1070               nlogchlsets = nlogchlsets + 1 
     1071 
     1072               CALL obs_rea_logchl( 0, logchldata(nlogchlsets), 1, & 
     1073                  &                 logchlfbfiles(jset:jset), & 
     1074                  &                 nlogchlvars, nlogchlextr, nitend-nit000+2, & 
     1075                  &                 dobsini, dobsend, ln_ignmis, .FALSE. ) 
     1076 
     1077               CALL obs_pre_logchl( logchldata(nlogchlsets), logchldatqc(nlogchlsets), & 
     1078                  &                 ln_logchl, ln_nea ) 
     1079             
     1080            ENDDO 
     1081 
     1082         ELSE                              ! Original file format 
     1083 
     1084            nlogchlsets = nlogchlsets + 1 
     1085 
     1086            CALL obs_rea_logchl( 1, logchldata(nlogchlsets), jnumlogchl, & 
     1087               &                 logchlfiles(1:jnumlogchl), & 
     1088               &                 nlogchlvars, nlogchlextr, nitend-nit000+2, & 
     1089               &                 dobsini, dobsend, ln_ignmis, .FALSE. ) 
     1090 
     1091            CALL obs_pre_logchl( logchldata(nlogchlsets), logchldatqc(nlogchlsets), & 
     1092               &                 ln_logchl, ln_nea ) 
     1093 
     1094         ENDIF 
     1095  
     1096      ENDIF 
    10011097      
    10021098   END SUBROUTINE dia_obs_init 
     
    10161112      !!               - Sea surface salinity 
    10171113      !!               - Velocity component (U,V) profiles 
     1114      !!               - Sea surface log10(chlorophyll) 
    10181115      !! 
    10191116      !! ** Action  :  
     
    10531150         & frld 
    10541151#endif 
     1152#if defined key_hadocc 
     1153      USE trc, ONLY :  &                ! HadOCC chlorophyll 
     1154         & HADOCC_CHL, & 
     1155         & HADOCC_FILL_FLT 
     1156#elif defined key_medusa && defined key_foam_medusa 
     1157      USE trc, ONLY :  &                ! MEDUSA chlorophyll 
     1158         & MEDUSA_CHL, & 
     1159         & MEDUSA_FILL_FLT 
     1160#elif defined key_fabm 
     1161      !USE ???                           ! ERSEM chlorophyll 
     1162#endif 
    10551163      IMPLICIT NONE 
    10561164 
     
    10641172      INTEGER :: jseaiceset             ! sea ice data set loop variable 
    10651173      INTEGER :: jveloset               ! velocity profile data loop variable 
     1174      INTEGER :: jlogchlset             ! logchl data set loop variable 
    10661175      INTEGER :: jvar                   ! Variable number     
    10671176#if ! defined key_lim2 && ! defined key_lim3 
    10681177      REAL(wp), POINTER, DIMENSION(:,:) :: frld    
    10691178#endif 
     1179      REAL(wp) :: tiny                  ! small number 
     1180      REAL(wp), DIMENSION(jpi,jpj) :: & 
     1181         logchl                         ! array for log chlorophyll 
     1182      REAL(wp), DIMENSION(jpi,jpj) :: & 
     1183         maskchl                        ! array for special chlorophyll mask 
    10701184      CHARACTER(LEN=20) :: datestr=" ",timestr=" " 
    10711185  
     
    11771291      ENDIF 
    11781292 
     1293      IF ( ln_logchl ) THEN 
     1294 
     1295#if defined key_hadocc 
     1296         logchl(:,:)  = HADOCC_CHL(:,:,1)    ! (not log) chlorophyll from HadOCC 
     1297#elif defined key_medusa && defined key_foam_medusa 
     1298         logchl(:,:)  = MEDUSA_CHL(:,:,1)    ! (not log) chlorophyll from HadOCC 
     1299#elif defined key_fabm 
     1300         !logchl(:,:)  =  ???                 ! (not log) chlorophyll from ERSEM 
     1301         CALL ctl_stop( ' Trying to run logchl observation operator', & 
     1302            &           ' but not properly implemented for FABM-ERSEM yet' ) 
     1303#else 
     1304         CALL ctl_stop( ' Trying to run logchl observation operator', & 
     1305            &           ' but no biogeochemical model appears to have been defined' ) 
     1306#endif 
     1307 
     1308         maskchl(:,:) = tmask(:,:,1)         ! create a special mask to exclude certain things 
     1309 
     1310         ! Take the log10 where we can, otherwise exclude 
     1311         tiny = 1.0e-20 
     1312         WHERE(logchl(:,:) > tiny .AND. logchl(:,:) /= obfillflt ) 
     1313            logchl(:,:)  = LOG10(logchl(:,:)) 
     1314         ELSEWHERE 
     1315            logchl(:,:)  = obfillflt 
     1316            maskchl(:,:) = 0 
     1317         END WHERE 
     1318 
     1319         DO jlogchlset = 1, nlogchlsets 
     1320             CALL obs_logchl_opt( logchldatqc(jlogchlset),             & 
     1321               &                  kstp, jpi, jpj, nit000, logchl(:,:), & 
     1322               &                  maskchl(:,:), n2dint ) 
     1323         END DO          
     1324      ENDIF  
     1325 
    11791326#if ! defined key_lim2 && ! defined key_lim3 
    11801327      CALL wrk_dealloc(jpi,jpj,frld)  
     
    12091356      INTEGER :: jsstset                  ! SST data set loop variable 
    12101357      INTEGER :: jseaiceset               ! Sea Ice data set loop variable 
     1358      INTEGER :: jlogchlset               ! logchl data set loop variable 
    12111359      INTEGER :: jset 
    12121360      INTEGER :: jfbini 
     
    14531601         ENDIF 
    14541602          
     1603      ENDIF 
     1604 
     1605      !  - log10(chlorophyll) 
     1606      IF ( ln_logchl ) THEN 
     1607 
     1608         ! Copy data from logchldatqc to logchldata structures 
     1609         DO jlogchlset = 1, nlogchlsets 
     1610 
     1611            CALL obs_surf_decompress( logchldatqc(jlogchlset), & 
     1612                 &                    logchldata(jlogchlset), .TRUE., numout ) 
     1613 
     1614         END DO 
     1615          
     1616         ! Mark as bad observations with no valid model counterpart due to activities in dia_obs 
     1617         ! Seem to need to set to fill value rather than marking as bad to be effective, so do both 
     1618         DO jlogchlset = 1, nlogchlsets 
     1619            WHERE ( logchldata(jlogchlset)%rmod(:,1) == obfillflt ) 
     1620               logchldata(jlogchlset)%nqc(:)    = 1 
     1621               logchldata(jlogchlset)%robs(:,1) = obfillflt 
     1622            END WHERE 
     1623         END DO 
     1624 
     1625         ! Write the logchl data 
     1626         DO jlogchlset = 1, nlogchlsets 
     1627       
     1628            WRITE(cdtmp,'(A,I2.2)')'logchlfb_',jlogchlset 
     1629            CALL obs_wri_logchl( cdtmp, logchldata(jlogchlset) ) 
     1630 
     1631         END DO 
     1632 
    14551633      ENDIF 
    14561634 
  • branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/obs_logchl.F90

    r6845 r6854  
    1 MODULE obs_seaice 
     1MODULE obs_logchl 
    22   !!===================================================================== 
    3    !!                       ***  MODULE  obs_seaice  *** 
    4    !! Observation diagnostics: Storage space for sea ice observations 
     3   !!                       ***  MODULE  obs_logchl  *** 
     4   !! Observation diagnostics: Storage space for logchl observations 
    55   !!                          arrays and additional flags etc. 
    66   !!===================================================================== 
     
    1313    
    1414   !! * Modules used  
    15    USE obs_surf_def ! Definition of sea ice data types and tools 
     15   USE obs_surf_def ! Definition of surface data types and tools 
    1616 
    1717   IMPLICIT NONE 
     
    2222   PRIVATE 
    2323 
    24    PUBLIC nseaicevars, nseaiceextr, nseaicesets, seaicedata, seaicedatqc 
     24   PUBLIC nlogchlvars, nlogchlextr, nlogchlsets, logchldata, logchldatqc 
    2525 
    2626   !! * Shared Module variables 
    27    INTEGER :: nseaicevars                               ! Number of seaicedata variables 
    28    INTEGER :: nseaiceextr                               ! Number of seaicedata extra  
    29                                                      ! variables 
    30    INTEGER :: nseaicesets                               ! Number of seaicedata sets 
    31    TYPE(obs_surf), POINTER, DIMENSION(:) :: seaicedata  ! Initial sea ice data 
    32    TYPE(obs_surf), POINTER, DIMENSION(:) :: seaicedatqc ! Sea ice data after quality control 
     27   INTEGER :: nlogchlvars                               ! Number of logchldata variables 
     28   INTEGER :: nlogchlextr                               ! Number of logchldata extra  
     29                                                        ! variables 
     30   INTEGER :: nlogchlsets                               ! Number of logchldata sets 
     31   TYPE(obs_surf), POINTER, DIMENSION(:) :: logchldata  ! Initial logchl data 
     32   TYPE(obs_surf), POINTER, DIMENSION(:) :: logchldatqc ! Sea ice data after quality control 
    3333 
    34 END MODULE obs_seaice 
     34END MODULE obs_logchl 
    3535 
  • branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/obs_logchl_io.F90

    r6845 r6854  
    1 MODULE obs_seaice_io 
     1MODULE obs_logchl_io 
    22   !!====================================================================== 
    3    !!                       ***  MODULE obs_seaice_io  *** 
    4    !! Observation operators : I/O for GHRSEAICE files 
     3   !!                       ***  MODULE obs_logchl_io  *** 
     4   !! Observation operators : I/O for logchl files 
    55   !!====================================================================== 
    66   !! History :  
     
    88   !!---------------------------------------------------------------------- 
    99   !!---------------------------------------------------------------------- 
    10    !!   read_seaicefile    :  Read a obfbdata structure from an GHRSEAICE file 
     10   !!   read_logchlfile    :  Read a obfbdata structure from a logchl file 
    1111   !!---------------------------------------------------------------------- 
    1212   USE par_kind 
     
    2626CONTAINS 
    2727 
    28 #include "obsseaice_io.h90" 
     28#include "obslogchl_io.h90" 
    2929 
    30 END MODULE obs_seaice_io 
     30END MODULE obs_logchl_io 
  • branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/obs_oper.F90

    r6301 r6854  
    2323   !!   obs_vel_opt :    Compute the model counterpart of zonal and meridional 
    2424   !!                    components of velocity from observations. 
     25   !!   obs_logchl_opt : Compute the model counterpart of log10(chlorophyll) 
     26   !!                    observations 
    2527   !!---------------------------------------------------------------------- 
    2628 
     
    6365      &   obs_sss_opt, &  ! Compute the model counterpart of SSS observations 
    6466      &   obs_seaice_opt, & 
    65       &   obs_vel_opt     ! Compute the model counterpart of velocity profile data 
     67      &   obs_vel_opt, &  ! Compute the model counterpart of velocity profile data 
     68      &   obs_logchl_opt  ! Compute the model counterpart of logchl data 
    6669 
    6770   INTEGER, PARAMETER, PUBLIC :: imaxavtypes = 20 ! Max number of daily avgd obs types 
     
    20522055   END SUBROUTINE obs_vel_opt 
    20532056 
     2057   SUBROUTINE obs_logchl_opt( logchldatqc, kt, kpi, kpj, kit000, & 
     2058      &                    plogchln, plogchlmask, k2dint ) 
     2059 
     2060      !!----------------------------------------------------------------------- 
     2061      !! 
     2062      !!                     ***  ROUTINE obs_logchl_opt  *** 
     2063      !! 
     2064      !! ** Purpose : Compute the model counterpart of logchl 
     2065      !!              data by interpolating from the model grid to the  
     2066      !!              observation point. 
     2067      !! 
     2068      !! ** Method  : Linearly interpolate to each observation point using  
     2069      !!              the model values at the corners of the surrounding grid box. 
     2070      !! 
     2071      !!    The now model logchl is first computed at the obs (lon, lat) point. 
     2072      !! 
     2073      !!    Several horizontal interpolation schemes are available: 
     2074      !!        - distance-weighted (great circle) (k2dint = 0) 
     2075      !!        - distance-weighted (small angle)  (k2dint = 1) 
     2076      !!        - bilinear (geographical grid)     (k2dint = 2) 
     2077      !!        - bilinear (quadrilateral grid)    (k2dint = 3) 
     2078      !!        - polynomial (quadrilateral grid)  (k2dint = 4) 
     2079      !! 
     2080      !! 
     2081      !! ** Action  : 
     2082      !! 
     2083      !! History : 
     2084      !!       
     2085      !!----------------------------------------------------------------------- 
     2086 
     2087      !! * Modules used 
     2088      USE obs_surf_def  ! Definition of storage space for surface observations 
     2089 
     2090      IMPLICIT NONE 
     2091 
     2092      !! * Arguments 
     2093      TYPE(obs_surf), INTENT(INOUT) :: logchldatqc     ! Subset of surface data not failing screening 
     2094      INTEGER, INTENT(IN) :: kt       ! Time step 
     2095      INTEGER, INTENT(IN) :: kpi      ! Model grid parameters 
     2096      INTEGER, INTENT(IN) :: kpj 
     2097      INTEGER, INTENT(IN) :: kit000   ! Number of the first time step  
     2098                                      !   (kit000-1 = restart time) 
     2099      INTEGER, INTENT(IN) :: k2dint   ! Horizontal interpolation type (see header) 
     2100      REAL(KIND=wp), INTENT(IN), DIMENSION(kpi,kpj) :: & 
     2101         & plogchln,  &    ! Model logchl field 
     2102         & plogchlmask     ! Land-sea mask 
     2103          
     2104      !! * Local declarations 
     2105      INTEGER :: ji 
     2106      INTEGER :: jj 
     2107      INTEGER :: jobs 
     2108      INTEGER :: inrc 
     2109      INTEGER :: ilogchl 
     2110      INTEGER :: iobs 
     2111        
     2112      REAL(KIND=wp) :: zlam 
     2113      REAL(KIND=wp) :: zphi 
     2114      REAL(KIND=wp) :: zext(1), zobsmask(1) 
     2115      REAL(kind=wp), DIMENSION(2,2,1) :: & 
     2116         & zweig 
     2117      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: & 
     2118         & zmask, & 
     2119         & zlogchll, & 
     2120         & zglam, & 
     2121         & zgphi 
     2122      INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: & 
     2123         & igrdi, & 
     2124         & igrdj 
     2125 
     2126      !------------------------------------------------------------------------ 
     2127      ! Local initialization  
     2128      !------------------------------------------------------------------------ 
     2129      ! ... Record and data counters 
     2130      inrc = kt - kit000 + 2 
     2131      ilogchl = logchldatqc%nsstp(inrc) 
     2132 
     2133      ! Get the data for interpolation 
     2134       
     2135      ALLOCATE( & 
     2136         & igrdi(2,2,ilogchl), & 
     2137         & igrdj(2,2,ilogchl), & 
     2138         & zglam(2,2,ilogchl), & 
     2139         & zgphi(2,2,ilogchl), & 
     2140         & zmask(2,2,ilogchl), & 
     2141         & zlogchll(2,2,ilogchl)  & 
     2142         & ) 
     2143       
     2144      DO jobs = logchldatqc%nsurfup + 1, logchldatqc%nsurfup + ilogchl 
     2145         iobs = jobs - logchldatqc%nsurfup 
     2146         igrdi(1,1,iobs) = logchldatqc%mi(jobs)-1 
     2147         igrdj(1,1,iobs) = logchldatqc%mj(jobs)-1 
     2148         igrdi(1,2,iobs) = logchldatqc%mi(jobs)-1 
     2149         igrdj(1,2,iobs) = logchldatqc%mj(jobs) 
     2150         igrdi(2,1,iobs) = logchldatqc%mi(jobs) 
     2151         igrdj(2,1,iobs) = logchldatqc%mj(jobs)-1 
     2152         igrdi(2,2,iobs) = logchldatqc%mi(jobs) 
     2153         igrdj(2,2,iobs) = logchldatqc%mj(jobs) 
     2154      END DO 
     2155       
     2156      CALL obs_int_comm_2d( 2, 2, ilogchl, & 
     2157         &                  igrdi, igrdj, glamt, zglam ) 
     2158      CALL obs_int_comm_2d( 2, 2, ilogchl, & 
     2159         &                  igrdi, igrdj, gphit, zgphi ) 
     2160      CALL obs_int_comm_2d( 2, 2, ilogchl, & 
     2161         &                  igrdi, igrdj, plogchlmask, zmask ) 
     2162      CALL obs_int_comm_2d( 2, 2, ilogchl, & 
     2163         &                  igrdi, igrdj, plogchln, zlogchll ) 
     2164       
     2165      DO jobs = logchldatqc%nsurfup + 1, logchldatqc%nsurfup + ilogchl 
     2166          
     2167         iobs = jobs - logchldatqc%nsurfup 
     2168          
     2169         IF ( kt /= logchldatqc%mstp(jobs) ) THEN 
     2170             
     2171            IF(lwp) THEN 
     2172               WRITE(numout,*) 
     2173               WRITE(numout,*) ' E R R O R : Observation',              & 
     2174                  &            ' time step is not consistent with the', & 
     2175                  &            ' model time step' 
     2176               WRITE(numout,*) ' =========' 
     2177               WRITE(numout,*) 
     2178               WRITE(numout,*) ' Record  = ', jobs,                & 
     2179                  &            ' kt      = ', kt,                  & 
     2180                  &            ' mstp    = ', logchldatqc%mstp(jobs), & 
     2181                  &            ' ntyp    = ', logchldatqc%ntyp(jobs) 
     2182            ENDIF 
     2183            CALL ctl_stop( 'obs_logchl_opt', 'Inconsistent time' ) 
     2184             
     2185         ENDIF 
     2186          
     2187         zlam = logchldatqc%rlam(jobs) 
     2188         zphi = logchldatqc%rphi(jobs) 
     2189          
     2190         ! Get weights to interpolate the model logchl to the observation point 
     2191         CALL obs_int_h2d_init( 1, 1, k2dint, zlam, zphi,         & 
     2192            &                   zglam(:,:,iobs), zgphi(:,:,iobs), & 
     2193            &                   zmask(:,:,iobs), zweig, zobsmask ) 
     2194          
     2195         ! ... Interpolate the model logchl to the observation point 
     2196         CALL obs_int_h2d( 1, 1,      & 
     2197            &              zweig, zlogchll(:,:,iobs),  zext ) 
     2198          
     2199         logchldatqc%rmod(jobs,1) = zext(1) 
     2200          
     2201      END DO 
     2202       
     2203      ! Deallocate the data for interpolation 
     2204      DEALLOCATE( & 
     2205         & igrdi,    & 
     2206         & igrdj,    & 
     2207         & zglam,    & 
     2208         & zgphi,    & 
     2209         & zmask,    & 
     2210         & zlogchll  & 
     2211         & ) 
     2212       
     2213      logchldatqc%nsurfup = logchldatqc%nsurfup + ilogchl 
     2214 
     2215   END SUBROUTINE obs_logchl_opt 
     2216 
    20542217END MODULE obs_oper 
    20552218 
  • branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/obs_prep.F90

    r6301 r6854  
    1212   !!   obs_pre_seaice : First level check and screening of sea ice observations 
    1313   !!   obs_pre_vel  : First level check and screening of velocity obs. 
     14   !!   obs_pre_logchl : First level check and screening of logchl obs. 
    1415   !!   obs_scr      : Basic screening of the observations 
    1516   !!   obs_coo_tim  : Compute number of time steps to the observation time 
     
    4142      & obs_pre_seaice, & ! First level check and screening of sea ice data 
    4243      & obs_pre_vel, &     ! First level check and screening of velocity profiles 
     44      & obs_pre_logchl, & ! First level check and screening of logchl data 
    4345      & calc_month_len     ! Calculate the number of days in the months of a year   
    4446 
     
    11861188   END SUBROUTINE obs_pre_vel 
    11871189 
     1190   SUBROUTINE obs_pre_logchl( logchldata, logchldatqc, ld_logchl, ld_nea ) 
     1191      !!---------------------------------------------------------------------- 
     1192      !!                    ***  ROUTINE obs_pre_logchl  *** 
     1193      !! 
     1194      !! ** Purpose : First level check and screening of logchl observations 
     1195      !! 
     1196      !! ** Method  : First level check and screening of logchl observations 
     1197      !! 
     1198      !! ** Action  :  
     1199      !! 
     1200      !! References : 
     1201      !!    
     1202      !! History : 
     1203      !!---------------------------------------------------------------------- 
     1204      !! * Modules used 
     1205      USE domstp              ! Domain: set the time-step 
     1206      USE par_oce             ! Ocean parameters 
     1207      USE dom_oce, ONLY : &   ! Geographical information 
     1208         & glamt,   & 
     1209         & gphit,   & 
     1210         & tmask 
     1211      !! * Arguments 
     1212      TYPE(obs_surf), INTENT(INOUT) :: logchldata     ! Full set of logchl data 
     1213      TYPE(obs_surf), INTENT(INOUT) :: logchldatqc    ! Subset of logchl data not failing screening 
     1214      LOGICAL :: ld_logchl     ! Switch for logchl data 
     1215      LOGICAL :: ld_nea        ! Switch for rejecting observation near land 
     1216      !! * Local declarations 
     1217      INTEGER :: iyea0         ! Initial date 
     1218      INTEGER :: imon0         !  - (year, month, day, hour, minute) 
     1219      INTEGER :: iday0     
     1220      INTEGER :: ihou0     
     1221      INTEGER :: imin0 
     1222      INTEGER :: icycle       ! Current assimilation cycle 
     1223                              ! Counters for observations that 
     1224      INTEGER :: iotdobs      !  - outside time domain 
     1225      INTEGER :: iosdsobs     !  - outside space domain 
     1226      INTEGER :: ilansobs     !  - within a model land cell 
     1227      INTEGER :: inlasobs     !  - close to land 
     1228      INTEGER :: igrdobs      !  - fail the grid search 
     1229                              ! Global counters for observations that 
     1230      INTEGER :: iotdobsmpp   !  - outside time domain 
     1231      INTEGER :: iosdsobsmpp  !  - outside space domain 
     1232      INTEGER :: ilansobsmpp  !  - within a model land cell 
     1233      INTEGER :: inlasobsmpp  !  - close to land 
     1234      INTEGER :: igrdobsmpp   !  - fail the grid search 
     1235      LOGICAL, DIMENSION(:), ALLOCATABLE :: &  
     1236         & llvalid            ! data selection 
     1237      INTEGER :: jobs         ! Obs. loop variable 
     1238      INTEGER :: jstp         ! Time loop variable 
     1239      INTEGER :: inrc         ! Time index variable 
     1240 
     1241      IF (lwp) WRITE(numout,*)'obs_pre_logchl : Preparing the logchl observations...' 
     1242 
     1243      ! Initial date initialization (year, month, day, hour, minute) 
     1244      iyea0 =   ndate0 / 10000 
     1245      imon0 = ( ndate0 - iyea0 * 10000 ) / 100 
     1246      iday0 =   ndate0 - iyea0 * 10000 - imon0 * 100 
     1247      ihou0 = 0 
     1248      imin0 = 0 
     1249 
     1250      icycle = no     ! Assimilation cycle 
     1251 
     1252      ! Diagnostics counters for various failures. 
     1253 
     1254      iotdobs  = 0 
     1255      igrdobs  = 0 
     1256      iosdsobs = 0 
     1257      ilansobs = 0 
     1258      inlasobs = 0 
     1259 
     1260      ! ----------------------------------------------------------------------- 
     1261      ! Find time coordinate for logchl data 
     1262      ! ----------------------------------------------------------------------- 
     1263 
     1264      CALL obs_coo_tim( icycle, & 
     1265         &              iyea0,   imon0,   iday0,   ihou0,   imin0,      & 
     1266         &              logchldata%nsurf,   logchldata%nyea, logchldata%nmon, & 
     1267         &              logchldata%nday,    logchldata%nhou, logchldata%nmin, & 
     1268         &              logchldata%nqc,     logchldata%mstp, iotdobs        ) 
     1269      CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp ) 
     1270      ! ----------------------------------------------------------------------- 
     1271      ! Check for logchl data failing the grid search 
     1272      ! ----------------------------------------------------------------------- 
     1273 
     1274      CALL obs_coo_grd( logchldata%nsurf,   logchldata%mi, logchldata%mj, & 
     1275         &              logchldata%nqc,     igrdobs                         ) 
     1276      CALL obs_mpp_sum_integer( igrdobs, igrdobsmpp ) 
     1277 
     1278      ! ----------------------------------------------------------------------- 
     1279      ! Check for land points.  
     1280      ! ----------------------------------------------------------------------- 
     1281 
     1282      CALL obs_coo_spc_2d( logchldata%nsurf,                 & 
     1283         &                 jpi,             jpj,             & 
     1284         &                 logchldata%mi,   logchldata%mj,   &  
     1285         &                 logchldata%rlam, logchldata%rphi, & 
     1286         &                 glamt,           gphit,           & 
     1287         &                 tmask(:,:,1),    logchldata%nqc,  & 
     1288         &                 iosdsobs,        ilansobs,        & 
     1289         &                 inlasobs,        ld_nea           )  
     1290          
     1291      CALL obs_mpp_sum_integer( iosdsobs, iosdsobsmpp ) 
     1292      CALL obs_mpp_sum_integer( ilansobs, ilansobsmpp ) 
     1293      CALL obs_mpp_sum_integer( inlasobs, inlasobsmpp ) 
     1294 
     1295      ! ----------------------------------------------------------------------- 
     1296      ! Copy useful data from the logchldata data structure to 
     1297      ! the logchldatqc data structure  
     1298      ! ----------------------------------------------------------------------- 
     1299 
     1300      ! Allocate the selection arrays 
     1301 
     1302      ALLOCATE( llvalid(logchldata%nsurf) ) 
     1303       
     1304      ! We want all data which has qc flags <= 0 
     1305 
     1306      llvalid(:)  = ( logchldata%nqc(:)  <= 10 ) 
     1307 
     1308      ! The actual copying 
     1309 
     1310      CALL obs_surf_compress( logchldata,     logchldatqc,       .TRUE.,  numout, & 
     1311         &                    lvalid=llvalid ) 
     1312 
     1313      ! Dellocate the selection arrays 
     1314      DEALLOCATE( llvalid ) 
     1315 
     1316      ! ----------------------------------------------------------------------- 
     1317      ! Print information about what observations are left after qc 
     1318      ! ----------------------------------------------------------------------- 
     1319 
     1320      ! Update the total observation counter array 
     1321       
     1322      IF(lwp) THEN 
     1323         WRITE(numout,*) 
     1324         WRITE(numout,*) 'obs_pre_logchl :' 
     1325         WRITE(numout,*) '~~~~~~~~~~~' 
     1326         WRITE(numout,*) 
     1327         WRITE(numout,*) ' logchl data outside time domain                  = ', & 
     1328            &            iotdobsmpp 
     1329         WRITE(numout,*) ' Remaining logchl data that failed grid search    = ', & 
     1330            &            igrdobsmpp 
     1331         WRITE(numout,*) ' Remaining logchl data outside space domain       = ', & 
     1332            &            iosdsobsmpp 
     1333         WRITE(numout,*) ' Remaining logchl data at land points             = ', & 
     1334            &            ilansobsmpp 
     1335         IF (ld_nea) THEN 
     1336            WRITE(numout,*) ' Remaining logchl data near land points (removed) = ', & 
     1337               &            inlasobsmpp 
     1338         ELSE 
     1339            WRITE(numout,*) ' Remaining logchl data near land points (kept)    = ', & 
     1340               &            inlasobsmpp 
     1341         ENDIF 
     1342         WRITE(numout,*) ' logchl data accepted                             = ', & 
     1343            &            logchldatqc%nsurfmpp 
     1344 
     1345         WRITE(numout,*) 
     1346         WRITE(numout,*) ' Number of observations per time step :' 
     1347         WRITE(numout,*) 
     1348         WRITE(numout,1997) 
     1349         WRITE(numout,1998) 
     1350      ENDIF 
     1351       
     1352      DO jobs = 1, logchldatqc%nsurf 
     1353         inrc = logchldatqc%mstp(jobs) + 2 - nit000 
     1354         logchldatqc%nsstp(inrc)  = logchldatqc%nsstp(inrc) + 1 
     1355      END DO 
     1356       
     1357      CALL obs_mpp_sum_integers( logchldatqc%nsstp, logchldatqc%nsstpmpp, & 
     1358         &                       nitend - nit000 + 2 ) 
     1359 
     1360      IF ( lwp ) THEN 
     1361         DO jstp = nit000 - 1, nitend 
     1362            inrc = jstp - nit000 + 2 
     1363            WRITE(numout,1999) jstp, logchldatqc%nsstpmpp(inrc) 
     1364         END DO 
     1365      ENDIF 
     1366 
     13671997  FORMAT(10X,'Time step',5X,'logchl data') 
     13681998  FORMAT(10X,'---------',5X,'------------') 
     13691999  FORMAT(10X,I9,5X,I17) 
     1370       
     1371   END SUBROUTINE obs_pre_logchl 
     1372 
    11881373   SUBROUTINE obs_coo_tim( kcycle, & 
    11891374      &                    kyea0,   kmon0,   kday0,   khou0,   kmin0,     & 
  • branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/obs_read_logchl.F90

    r6845 r6854  
    1 MODULE obs_read_seaice 
     1MODULE obs_read_logchl 
    22   !!====================================================================== 
    3    !!                       ***  MODULE obs_read_seaice  *** 
    4    !! Observation diagnostics: Read the along track SEAICE data from 
    5    !!                          GHRSST or any SEAICE data from feedback files 
     3   !!                       ***  MODULE obs_read_logchl  *** 
     4   !! Observation diagnostics: Read the along track logchl data from 
     5   !!                          GHRSST or any logchl data from feedback files 
    66   !!====================================================================== 
    77 
    88   !!---------------------------------------------------------------------- 
    9    !!   obs_rea_seaice : Driver for reading seaice data from the GHRSST/feedback 
     9   !!   obs_rea_logchl : Driver for reading logchl data from the feedback 
    1010   !!---------------------------------------------------------------------- 
    1111 
     
    2121   USE obs_surf_def             ! Surface observation definitions 
    2222   USE obs_types                ! Observation type definitions 
    23    USE obs_seaice_io            ! I/O for seaice files 
    24    USE iom                      ! I/O of fields for Reynolds data 
     23   USE obs_logchl_io            ! I/O for logchl files 
     24   USE iom                      ! I/O 
    2525   USE netcdf                   ! NetCDF library 
    2626 
     
    3030   PRIVATE 
    3131 
    32    PUBLIC obs_rea_seaice      ! Read the seaice observations from the point data 
     32   PUBLIC obs_rea_logchl      ! Read the logchl observations from the point data 
    3333    
    3434   !!---------------------------------------------------------------------- 
     
    4040CONTAINS 
    4141 
    42    SUBROUTINE obs_rea_seaice( kformat, & 
    43       &                    seaicedata, knumfiles, cfilenames, & 
     42   SUBROUTINE obs_rea_logchl( kformat, & 
     43      &                    logchldata, knumfiles, cfilenames, & 
    4444      &                    kvars, kextr, kstp, ddobsini, ddobsend, & 
    4545      &                    ldignmis, ldmod ) 
    4646      !!--------------------------------------------------------------------- 
    4747      !! 
    48       !!                   *** ROUTINE obs_rea_seaice *** 
    49       !! 
    50       !! ** Purpose : Read from file the seaice data 
    51       !! 
    52       !! ** Method  : Depending on kformat either AVISO or 
     48      !!                   *** ROUTINE obs_rea_logchl *** 
     49      !! 
     50      !! ** Purpose : Read from file the logchl data 
     51      !! 
     52      !! ** Method  : Depending on kformat either old or new style 
    5353      !!              feedback data files are read 
    5454      !! 
     
    6363      !! * Arguments 
    6464      INTEGER :: kformat   ! Format of input data  
    65       !                    ! 0: Feedback 
    66       !                    ! 1: GHRSST 
     65      !                    ! 0: New-style feedback 
     66      !                    ! 1: Old-style feedback 
    6767      TYPE(obs_surf), INTENT(INOUT) :: & 
    68          & seaicedata     ! seaice data to be read 
     68         & logchldata     ! logchl data to be read 
    6969      INTEGER, INTENT(IN) :: knumfiles   ! Number of corio format files to read in 
    7070      CHARACTER(LEN=128), INTENT(IN) :: cfilenames(knumfiles) ! File names to read in 
    71       INTEGER, INTENT(IN) :: kvars    ! Number of variables in seaicedata 
    72       INTEGER, INTENT(IN) :: kextr    ! Number of extra fields for each var in seaicedata 
     71      INTEGER, INTENT(IN) :: kvars    ! Number of variables in logchldata 
     72      INTEGER, INTENT(IN) :: kextr    ! Number of extra fields for each var in logchldata 
    7373      INTEGER, INTENT(IN) :: kstp     ! Ocean time-step index 
    7474      LOGICAL, INTENT(IN) :: ldignmis   ! Ignore missing files 
     
    7878          
    7979      !! * Local declarations 
    80       CHARACTER(LEN=14), PARAMETER :: cpname='obs_rea_seaice' 
     80      CHARACTER(LEN=14), PARAMETER :: cpname='obs_rea_logchl' 
    8181      INTEGER :: ji 
    8282      INTEGER :: jj 
     
    9595         & irefdate 
    9696      INTEGER :: iobsmpp 
    97       INTEGER, PARAMETER :: iseaicemaxtype = 1024 
    98       INTEGER, DIMENSION(0:iseaicemaxtype) :: & 
     97      INTEGER, PARAMETER :: ilogchlmaxtype = 1024 
     98      INTEGER, DIMENSION(0:ilogchlmaxtype) :: & 
    9999         & ityp, & 
    100100         & itypmpp 
     
    105105         & iindx,    & 
    106106         & ifileidx, & 
    107          & iseaiceidx 
     107         & ilogchlidx 
    108108      INTEGER :: itype 
    109109      REAL(wp), DIMENSION(:), ALLOCATABLE :: & 
     
    143143      ALLOCATE( inpfiles(inobf) ) 
    144144 
    145       seaice_files : DO jj = 1, inobf 
     145      logchl_files : DO jj = 1, inobf 
    146146           
    147147         !--------------------------------------------------------------------- 
     
    150150         IF(lwp) THEN 
    151151            WRITE(numout,*) 
    152             WRITE(numout,*) ' obs_rea_seaice : Reading from file = ', & 
     152            WRITE(numout,*) ' obs_rea_logchl : Reading from file = ', & 
    153153               & TRIM( TRIM( cfilenames(jj) ) ) 
    154154            WRITE(numout,*) ' ~~~~~~~~~~~~~~' 
     
    177177             
    178178            !------------------------------------------------------------------ 
    179             !  Close the file since it is opened in read_proffile 
     179            !  Close the file since it is opened elsewhere 
    180180            !------------------------------------------------------------------ 
    181181             
     
    183183 
    184184            !------------------------------------------------------------------ 
    185             !  Read the profile file into inpfiles 
     185            !  Read the file into inpfiles 
    186186            !------------------------------------------------------------------ 
    187187            IF ( kformat == 0 ) THEN 
     
    200200               ENDIF 
    201201            ELSEIF ( kformat == 1) THEN 
    202                CALL read_seaice( TRIM( cfilenames(jj) ), inpfiles(jj), & 
     202               CALL read_logchl( TRIM( cfilenames(jj) ), inpfiles(jj), & 
    203203               &                 numout, lwp, .TRUE. ) 
    204204            ELSE 
     
    291291         ENDIF 
    292292 
    293       END DO seaice_files 
     293      END DO logchl_files 
    294294 
    295295      !----------------------------------------------------------------------- 
     
    298298 
    299299      !--------------------------------------------------------------------- 
    300       !  Loop over input data files to count total number of profiles 
     300      !  Loop over input data files to count total number of obs 
    301301      !--------------------------------------------------------------------- 
    302302      iobstot = 0 
     
    311311 
    312312      ALLOCATE( iindx(iobstot), ifileidx(iobstot), & 
    313          &      iseaiceidx(iobstot), zdat(iobstot) ) 
     313         &      ilogchlidx(iobstot), zdat(iobstot) ) 
    314314      jk = 0 
    315315      DO jj = 1, inobf 
     
    319319               jk = jk + 1 
    320320               ifileidx(jk) = jj 
    321                iseaiceidx(jk) = ji 
     321               ilogchlidx(jk) = ji 
    322322               zdat(jk)     = inpfiles(jj)%ptim(ji) 
    323323            ENDIF 
     
    328328         &               iindx   ) 
    329329       
    330       CALL obs_surf_alloc( seaicedata, iobs, &  
     330      CALL obs_surf_alloc( logchldata, iobs, &  
    331331                           kvars, kextr, kstp, jpi, jpj ) 
    332332       
    333       ! * Read obs/positions, QC, all variable and assign to seaicedata 
     333      ! * Read obs/positions, QC, all variable and assign to logchldata 
    334334  
    335335      iobs = 0 
     
    343343          
    344344         jj = ifileidx(iindx(jk)) 
    345          ji = iseaiceidx(iindx(jk)) 
     345         ji = ilogchlidx(iindx(jk)) 
    346346         IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND.  & 
    347347            & ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN 
     
    370370 
    371371 
    372                ! seaice time coordinates 
    373                seaicedata%nyea(iobs) = iyea 
    374                seaicedata%nmon(iobs) = imon 
    375                seaicedata%nday(iobs) = iday 
    376                seaicedata%nhou(iobs) = ihou 
    377                seaicedata%nmin(iobs) = imin 
     372               ! logchl time coordinates 
     373               logchldata%nyea(iobs) = iyea 
     374               logchldata%nmon(iobs) = imon 
     375               logchldata%nday(iobs) = iday 
     376               logchldata%nhou(iobs) = ihou 
     377               logchldata%nmin(iobs) = imin 
    378378                
    379                ! seaice space coordinates 
    380                seaicedata%rlam(iobs) = inpfiles(jj)%plam(ji) 
    381                seaicedata%rphi(iobs) = inpfiles(jj)%pphi(ji) 
     379               ! logchl space coordinates 
     380               logchldata%rlam(iobs) = inpfiles(jj)%plam(ji) 
     381               logchldata%rphi(iobs) = inpfiles(jj)%pphi(ji) 
    382382 
    383383               ! Coordinate search parameters 
    384                seaicedata%mi  (iobs) = inpfiles(jj)%iobsi(ji,1) 
    385                seaicedata%mj  (iobs) = inpfiles(jj)%iobsj(ji,1) 
     384               logchldata%mi  (iobs) = inpfiles(jj)%iobsi(ji,1) 
     385               logchldata%mj  (iobs) = inpfiles(jj)%iobsj(ji,1) 
    386386                
    387387               ! Instrument type 
     
    392392                  itype = 0 
    393393               ENDIF 
    394                seaicedata%ntyp(iobs) = itype 
    395                IF ( itype < iseaicemaxtype + 1 ) THEN 
     394               logchldata%ntyp(iobs) = itype 
     395               IF ( itype < ilogchlmaxtype + 1 ) THEN 
    396396                  ityp(itype+1) = ityp(itype+1) + 1 
    397397               ELSE 
    398                   IF(lwp)WRITE(numout,*)'WARNING:Increase iseaicemaxtype in ',& 
     398                  IF(lwp)WRITE(numout,*)'WARNING:Increase ilogchlmaxtype in ',& 
    399399                     &                  cpname 
    400400               ENDIF 
    401401 
    402402               ! Bookkeeping data to match observations 
    403                seaicedata%nsidx(iobs) = iobs 
    404                seaicedata%nsfil(iobs) = iindx(jk) 
     403               logchldata%nsidx(iobs) = iobs 
     404               logchldata%nsfil(iobs) = iindx(jk) 
    405405 
    406406               ! QC flags 
    407                seaicedata%nqc(iobs) = inpfiles(jj)%ivqc(ji,1) 
     407               logchldata%nqc(iobs) = inpfiles(jj)%ivqc(ji,1) 
    408408 
    409409               ! Observed value 
    410                seaicedata%robs(iobs,1) = inpfiles(jj)%pob(1,ji,1) 
     410               logchldata%robs(iobs,1) = inpfiles(jj)%pob(1,ji,1) 
    411411 
    412412 
    413413               ! Model and MDT is set to fbrmdi unless read from file 
    414414               IF ( ldmod ) THEN 
    415                   seaicedata%rmod(iobs,1) = inpfiles(jj)%padd(1,ji,1,1) 
     415                  logchldata%rmod(iobs,1) = inpfiles(jj)%padd(1,ji,1,1) 
    416416               ELSE 
    417                   seaicedata%rmod(iobs,1) = fbrmdi 
     417                  logchldata%rmod(iobs,1) = fbrmdi 
    418418               ENDIF 
    419419            ENDIF 
     
    434434 
    435435         WRITE(numout,*) 
    436          WRITE(numout,'(1X,A)')'Seaice data types' 
     436         WRITE(numout,'(1X,A)')'logchl data types' 
    437437         WRITE(numout,'(1X,A)')'-----------------' 
    438438         DO jj = 1,8 
     
    450450      ! Deallocate temporary data 
    451451      !----------------------------------------------------------------------- 
    452       DEALLOCATE( ifileidx, iseaiceidx, zdat ) 
     452      DEALLOCATE( ifileidx, ilogchlidx, zdat ) 
    453453 
    454454      !----------------------------------------------------------------------- 
     
    460460      DEALLOCATE( inpfiles ) 
    461461 
    462    END SUBROUTINE obs_rea_seaice 
    463  
    464 END MODULE obs_read_seaice 
    465  
     462   END SUBROUTINE obs_rea_logchl 
     463 
     464END MODULE obs_read_logchl 
     465 
  • branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/obs_write.F90

    r5838 r6854  
    1111   !!   obs_wri_seaice: Write seaice observation related diagnostics 
    1212   !!   obs_wri_vel   : Write velocity observation diagnostics in NetCDF format 
     13   !!   obs_wri_logchl: Write logchl observation related diagnostics 
    1314   !!   obs_wri_stats : Print basic statistics on the data being written out 
    1415   !!---------------------------------------------------------------------- 
     
    4546      &   obs_wri_seaice, & ! Write seaice observation related diagnostics 
    4647      &   obs_wri_vel, &    ! Write velocity observation related diagnostics 
     48      &   obs_wri_logchl, & ! Write logchl observation related diagnostics 
    4749      &   obswriinfo 
    4850    
     
    930932   END SUBROUTINE obs_wri_vel 
    931933 
     934   SUBROUTINE obs_wri_logchl( cprefix, logchldata, padd, pext ) 
     935      !!----------------------------------------------------------------------- 
     936      !! 
     937      !!                     *** ROUTINE obs_wri_logchl  *** 
     938      !! 
     939      !! ** Purpose : Write logchl observation diagnostics 
     940      !!              related  
     941      !! 
     942      !! ** Method  : NetCDF 
     943      !!  
     944      !! ** Action  : 
     945      !! 
     946      !!----------------------------------------------------------------------- 
     947 
     948      !! * Modules used 
     949      IMPLICIT NONE 
     950 
     951      !! * Arguments 
     952      CHARACTER(LEN=*), INTENT(IN) :: cprefix       ! Prefix for output files 
     953      TYPE(obs_surf), INTENT(INOUT) :: logchldata   ! Full set of logchl 
     954      TYPE(obswriinfo), OPTIONAL :: padd            ! Additional info for each variable 
     955      TYPE(obswriinfo), OPTIONAL :: pext            ! Extra info 
     956 
     957      !! * Local declarations  
     958      TYPE(obfbdata) :: fbdata 
     959      CHARACTER(LEN=40) :: cfname             ! netCDF filename 
     960      CHARACTER(LEN=12), PARAMETER :: cpname = 'obs_wri_logchl' 
     961      INTEGER :: jo 
     962      INTEGER :: ja 
     963      INTEGER :: je 
     964      INTEGER :: nadd 
     965      INTEGER :: next 
     966 
     967      IF ( PRESENT( padd ) ) THEN 
     968         nadd = padd%inum 
     969      ELSE 
     970         nadd = 0 
     971      ENDIF 
     972 
     973      IF ( PRESENT( pext ) ) THEN 
     974         next = pext%inum 
     975      ELSE 
     976         next = 0 
     977      ENDIF 
     978 
     979      CALL init_obfbdata( fbdata ) 
     980 
     981      CALL alloc_obfbdata( fbdata, 1, logchldata%nsurf, 1, & 
     982         &                 1 + nadd, next, .TRUE. ) 
     983 
     984      fbdata%cname(1)      = 'LOGCHL' 
     985      fbdata%coblong(1)    = 'logchl concentration' 
     986      fbdata%cobunit(1)    = 'mg/m3' 
     987      DO je = 1, next 
     988         fbdata%cextname(je) = pext%cdname(je) 
     989         fbdata%cextlong(je) = pext%cdlong(je,1) 
     990         fbdata%cextunit(je) = pext%cdunit(je,1) 
     991      END DO 
     992      fbdata%caddname(1)   = 'Hx' 
     993      fbdata%caddlong(1,1) = 'Model interpolated LOGCHL' 
     994      fbdata%caddunit(1,1) = 'mg/m3' 
     995      fbdata%cgrid(1)      = 'T' 
     996      DO ja = 1, nadd 
     997         fbdata%caddname(1+ja) = padd%cdname(ja) 
     998         fbdata%caddlong(1+ja,1) = padd%cdlong(ja,1) 
     999         fbdata%caddunit(1+ja,1) = padd%cdunit(ja,1) 
     1000      END DO 
     1001 
     1002      WRITE(cfname, FMT="(A,'_fdbk_',I4.4,'.nc')") TRIM(cprefix), nproc 
     1003 
     1004      IF(lwp) THEN 
     1005         WRITE(numout,*) 
     1006         WRITE(numout,*)'obs_wri_logchl :' 
     1007         WRITE(numout,*)'~~~~~~~~~~~~~~~~' 
     1008         WRITE(numout,*)'Writing logchl feedback file : ',TRIM(cfname) 
     1009      ENDIF 
     1010 
     1011      ! Transform obs_prof data structure into obfbdata structure 
     1012      fbdata%cdjuldref = '19500101000000' 
     1013      DO jo = 1, logchldata%nsurf 
     1014         fbdata%plam(jo)      = logchldata%rlam(jo) 
     1015         fbdata%pphi(jo)      = logchldata%rphi(jo) 
     1016         WRITE(fbdata%cdtyp(jo),'(I4)') logchldata%ntyp(jo) 
     1017         fbdata%ivqc(jo,:)    = 0 
     1018         fbdata%ivqcf(:,jo,:) = 0 
     1019         IF ( logchldata%nqc(jo) > 10 ) THEN 
     1020            fbdata%ioqc(jo)    = 4 
     1021            fbdata%ioqcf(1,jo) = 0 
     1022            fbdata%ioqcf(2,jo) = logchldata%nqc(jo) - 10 
     1023         ELSE 
     1024            fbdata%ioqc(jo)    = MAX(logchldata%nqc(jo),1) 
     1025            fbdata%ioqcf(:,jo) = 0 
     1026         ENDIF 
     1027         fbdata%ipqc(jo)      = 0 
     1028         fbdata%ipqcf(:,jo)   = 0 
     1029         fbdata%itqc(jo)      = 0 
     1030         fbdata%itqcf(:,jo)   = 0 
     1031         fbdata%cdwmo(jo)     = '' 
     1032         fbdata%kindex(jo)    = logchldata%nsfil(jo) 
     1033         IF (ln_grid_global) THEN 
     1034            fbdata%iobsi(jo,1) = logchldata%mi(jo) 
     1035            fbdata%iobsj(jo,1) = logchldata%mj(jo) 
     1036         ELSE 
     1037            fbdata%iobsi(jo,1) = mig(logchldata%mi(jo)) 
     1038            fbdata%iobsj(jo,1) = mjg(logchldata%mj(jo)) 
     1039         ENDIF 
     1040         CALL greg2jul( 0, & 
     1041            &           logchldata%nmin(jo), & 
     1042            &           logchldata%nhou(jo), & 
     1043            &           logchldata%nday(jo), & 
     1044            &           logchldata%nmon(jo), & 
     1045            &           logchldata%nyea(jo), & 
     1046            &           fbdata%ptim(jo),   & 
     1047            &           krefdate = 19500101 ) 
     1048         fbdata%padd(1,jo,1,1) = logchldata%rmod(jo,1) 
     1049         fbdata%pob(1,jo,1)    = logchldata%robs(jo,1) 
     1050         fbdata%pdep(1,jo)     = 0.0 
     1051         fbdata%idqc(1,jo)     = 0 
     1052         fbdata%idqcf(:,1,jo)  = 0 
     1053         IF ( logchldata%nqc(jo) > 10 ) THEN 
     1054            fbdata%ivlqc(1,jo,1) = 4 
     1055            fbdata%ivlqcf(1,1,jo,1) = 0 
     1056            fbdata%ivlqcf(2,1,jo,1) = logchldata%nqc(jo) - 10 
     1057         ELSE 
     1058            fbdata%ivlqc(1,jo,1) = MAX(logchldata%nqc(jo),1) 
     1059            fbdata%ivlqcf(:,1,jo,1) = 0 
     1060         ENDIF 
     1061         fbdata%iobsk(1,jo,1)  = 0 
     1062         DO ja = 1, nadd 
     1063            fbdata%padd(1,jo,1+ja,1) = & 
     1064               & logchldata%rext(jo,padd%ipoint(ja)) 
     1065         END DO 
     1066         DO je = 1, next 
     1067            fbdata%pext(1,jo,je) = & 
     1068               & logchldata%rext(jo,pext%ipoint(je)) 
     1069         END DO 
     1070 
     1071      END DO 
     1072 
     1073      ! Write the obfbdata structure 
     1074      CALL write_obfbdata( cfname, fbdata ) 
     1075       
     1076      ! Output some basic statistics 
     1077      CALL obs_wri_stats( fbdata ) 
     1078 
     1079      CALL dealloc_obfbdata( fbdata ) 
     1080 
     1081   END SUBROUTINE obs_wri_logchl 
     1082 
    9321083   SUBROUTINE obs_wri_stats( fbdata ) 
    9331084      !!----------------------------------------------------------------------- 
  • branches/UKMO/dev_r4650_general_vert_coord_obsoper_surf_bgc/NEMOGCM/NEMO/OPA_SRC/OBS/obslogchl_io.h90

    r6845 r6854  
    55   !!---------------------------------------------------------------------- 
    66 
    7    SUBROUTINE read_seaice( cdfilename, inpfile, kunit, ldwp, ldgrid ) 
     7   SUBROUTINE read_logchl( cdfilename, inpfile, kunit, ldwp, ldgrid ) 
    88      !!--------------------------------------------------------------------- 
    99      !! 
    10       !!                     ** ROUTINE read_seaice ** 
    11       !! 
    12       !! ** Purpose : Read from file the SEAICE observations. 
     10      !!                     ** ROUTINE read_logchl ** 
     11      !! 
     12      !! ** Purpose : Read from file the logchl observations. 
    1313      !! 
    1414      !! ** Method  : The data file is a NetCDF file.  
     
    2828      LOGICAL          :: ldgrid     ! Save grid info in data structure 
    2929      !! * Local declarations 
    30       CHARACTER(LEN=12),PARAMETER :: cpname = 'read_seaice' 
     30      CHARACTER(LEN=12),PARAMETER :: cpname = 'read_logchl' 
    3131      INTEGER :: i_file_id    ! netcdf IDS 
    3232      INTEGER :: i_time_id 
     
    4141         & i_dtime, &       ! Offset in seconds since reference time 
    4242         & i_qc,    &       ! Quality control flag. 
    43          & i_type           ! Type of seaice measurement.             
     43         & i_type           ! Type of logchl measurement.             
    4444      REAL(wp), DIMENSION(:), POINTER :: & 
    4545         & z_phi,   &       ! Latitudes 
    4646         & z_lam            ! Longitudes 
    4747      REAL(wp), DIMENSION(:,:), POINTER :: & 
    48          & z_seaice         ! Seaice data      
     48         & z_logchl         ! logchl data      
    4949      INTEGER, PARAMETER :: imaxdim = 2    ! Assumed maximum for no. dims. in file 
    5050      INTEGER, DIMENSION(2) :: idims       ! Dimensions in file 
     
    9494         & z_phi        ( i_data                 ), &    
    9595         & z_lam        ( i_data                 ), &   
    96          & z_seaice        ( i_data,i_time  )  & 
     96         & z_logchl        ( i_data,i_time  )  & 
    9797         & ) 
    9898       
     
    124124      ! Get list of times for each ob in seconds relative to reference time 
    125125       
    126       CALL chkerr( nf90_inq_varid( i_file_id, 'SeaIce_dtime', i_var_id ), &  
     126      CALL chkerr( nf90_inq_varid( i_file_id, 'LogChl_dtime', i_var_id ), &  
    127127         &       cpname, __LINE__ ) 
    128128      idims(1) = i_data 
     
    164164         &         cpname, __LINE__ ) 
    165165       
    166       ! Get seaice data 
    167        
    168       CALL chkerr( nf90_inq_varid( i_file_id, 'sea_ice_concentration', & 
     166      ! Get logchl data 
     167       
     168      CALL chkerr( nf90_inq_varid( i_file_id, 'LogChl', & 
    169169         &                         i_var_id ), &  
    170170         &         cpname, __LINE__ ) 
     
    172172      idims(2) = i_time 
    173173      CALL chkdim( i_file_id, i_var_id, 2, idims, cpname, __LINE__ ) 
    174       CALL chkerr( nf90_get_var  ( i_file_id, i_var_id, z_seaice), & 
     174      CALL chkerr( nf90_get_var  ( i_file_id, i_var_id, z_logchl), & 
    175175         &       cpname, __LINE__ ) 
    176176      zoff = 0. 
    177       IF (nf90_inquire_attribute( i_file_id, i_var_id, "scale_factor") & 
    178          &                      == nf90_noerr) THEN 
    179          CALL chkerr( nf90_get_att( i_file_id, i_var_id, & 
    180             &                     "scale_factor",zsca), cpname, __LINE__ ) 
     177      IF (nf90_inquire_attribute( i_file_id, i_var_id, "add_offset") & 
     178         &                      == nf90_noerr) THEN 
     179         CALL chkerr( nf90_get_att( i_file_id, i_var_id, & 
     180            &                     "add_offset", zoff), cpname, __LINE__ ) 
    181181      ENDIF 
    182182      zsca = 1.0 
     
    192192            &                       "_FillValue",zfill), cpname, __LINE__ ) 
    193193      ENDIF 
    194       WHERE(z_seaice(:,:) /=  zfill) 
    195          z_seaice(:,:) = (zsca * z_seaice(:,:)) + zoff 
     194      WHERE(z_logchl(:,:) /=  zfill) 
     195         z_logchl(:,:) = (zsca * z_logchl(:,:)) + zoff 
    196196      ELSEWHERE 
    197          z_seaice(:,:) = fbrmdi 
     197         z_logchl(:,:) = fbrmdi 
    198198      END WHERE 
    199199       
     
    208208            &       cpname, __LINE__ ) 
    209209       
    210       ! Get seaice obs type 
     210      ! Get logchl obs type 
    211211       
    212212      i_type(:,:)=1 
     
    223223      CALL init_obfbdata( inpfile ) 
    224224      CALL alloc_obfbdata( inpfile, 1, iobs, 1, 0, 0, ldgrid ) 
    225       inpfile%cname(1) = 'SEAICE' 
     225      inpfile%cname(1) = 'LOGCHL' 
    226226 
    227227      ! Fill the obfbdata structure from input data 
     
    233233            iobs = iobs + 1 
    234234            ! Characters 
    235             WRITE(inpfile%cdwmo(iobs),'(A6,A2)') 'seaice','  ' 
     235            WRITE(inpfile%cdwmo(iobs),'(A6,A2)') 'logchl','  ' 
    236236            WRITE(inpfile%cdtyp(iobs),'(I4)') i_type(jobs,jtim) 
    237237            ! Real values 
    238238            inpfile%plam(iobs)         = z_lam(jobs) 
    239239            inpfile%pphi(iobs)         = z_phi(jobs) 
    240             inpfile%pob(1,iobs,1)      = z_seaice(jobs,jtim) 
     240            inpfile%pob(1,iobs,1)      = z_logchl(jobs,jtim) 
    241241            inpfile%ptim(iobs)         = & 
    242242               & REAL(i_reftime(jtim))/(60.*60.*24.) + & 
     
    245245            ! Integers 
    246246            inpfile%kindex(iobs)       = iobs 
    247             IF ( z_seaice(jobs,jtim) == fbrmdi ) THEN 
     247            IF ( z_logchl(jobs,jtim) == fbrmdi ) THEN 
    248248               inpfile%ioqc(iobs)      = 4 
    249249               inpfile%ivqc(iobs,1)    = 4  
     
    266266      END DO 
    267267 
    268    END SUBROUTINE read_seaice 
    269  
    270  
     268   END SUBROUTINE read_logchl 
     269 
     270 
Note: See TracChangeset for help on using the changeset viewer.