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 5998 for branches/2015/dev_MetOffice_merge_2015/NEMOGCM/NEMO/OPA_SRC/OBS/obs_read_prof.F90 – NEMO

Ignore:
Timestamp:
2015-12-04T11:56:46+01:00 (8 years ago)
Author:
timgraham
Message:

Merge in changes from OBS simplification branch (branches/2015/dev_r5072_UKMO2_OBS_simplification)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_MetOffice_merge_2015/NEMOGCM/NEMO/OPA_SRC/OBS/obs_read_prof.F90

    r4990 r5998  
    2525   USE netcdf                   ! NetCDF library 
    2626   USE obs_oper                 ! Observation operators 
    27    USE obs_prof_io              ! Profile files I/O (non-FB files) 
    2827   USE lib_mpp                  ! For ctl_warn/stop 
     28   USE obs_fbm                  ! Feedback routines 
    2929 
    3030   IMPLICIT NONE 
     
    3333   PRIVATE 
    3434 
    35    PUBLIC obs_rea_pro_dri  ! Read the profile observations  
     35   PUBLIC obs_rea_prof  ! Read the profile observations  
    3636 
    3737   !!---------------------------------------------------------------------- 
     
    4242 
    4343CONTAINS 
    44   
    45    SUBROUTINE obs_rea_pro_dri( kformat, & 
    46       &                        profdata, knumfiles, cfilenames, & 
    47       &                        kvars, kextr, kstp, ddobsini, ddobsend, & 
    48       &                        ldt3d, lds3d, ldignmis, ldsatt, ldavtimset, & 
    49       &                        ldmod, kdailyavtypes ) 
     44 
     45   SUBROUTINE obs_rea_prof( profdata, knumfiles, cdfilenames, & 
     46      &                     kvars, kextr, kstp, ddobsini, ddobsend, & 
     47      &                     ldvar1, ldvar2, ldignmis, ldsatt, & 
     48      &                     ldmod, kdailyavtypes ) 
    5049      !!--------------------------------------------------------------------- 
    5150      !! 
    52       !!                   *** ROUTINE obs_rea_pro_dri *** 
     51      !!                   *** ROUTINE obs_rea_prof *** 
    5352      !! 
    5453      !! ** Purpose : Read from file the profile observations 
    5554      !! 
    56       !! ** Method  : Depending on kformat either ENACT, CORIOLIS or 
    57       !!              feedback data files are read 
     55      !! ** Method  : Read feedback data in and transform to NEMO internal  
     56      !!              profile data structure 
    5857      !! 
    5958      !! ** Action  :  
     
    6362      !! History :   
    6463      !!      ! :  2009-09 (K. Mogensen) : New merged version of old routines 
     64      !!      ! :  2015-08 (M. Martin) : Merged profile and velocity routines 
    6565      !!---------------------------------------------------------------------- 
    66       !! * Modules used 
    67     
     66 
    6867      !! * Arguments 
    69       INTEGER ::  kformat    ! Format of input data 
    70       !                      ! 1: ENACT 
    71       !                      ! 2: Coriolis 
    72       TYPE(obs_prof), INTENT(OUT) ::  profdata     ! Profile data to be read 
    73       INTEGER, INTENT(IN) :: knumfiles      ! Number of files to read in 
     68      TYPE(obs_prof), INTENT(OUT) :: & 
     69         & profdata                     ! Profile data to be read 
     70      INTEGER, INTENT(IN) :: knumfiles  ! Number of files to read 
    7471      CHARACTER(LEN=128), INTENT(IN) ::  & 
    75          & cfilenames(knumfiles)  ! File names to read in 
     72         & cdfilenames(knumfiles)        ! File names to read in 
    7673      INTEGER, INTENT(IN) :: kvars      ! Number of variables in profdata 
    77       INTEGER, INTENT(IN) :: kextr      ! Number of extra fields for each var in profdata 
    78       INTEGER, INTENT(IN) :: kstp        ! Ocean time-step index 
    79       LOGICAL, INTENT(IN) :: ldt3d       ! Observed variables switches 
    80       LOGICAL, INTENT(IN) :: lds3d 
    81       LOGICAL, INTENT(IN) :: ldignmis    ! Ignore missing files 
    82       LOGICAL, INTENT(IN) :: ldsatt      ! Compute salinity at all temperature points 
    83       LOGICAL, INTENT(IN) :: ldavtimset  ! Correct time for daily averaged data 
    84       LOGICAL, INTENT(IN) :: ldmod       ! Initialize model from input data 
    85       REAL(KIND=dp), INTENT(IN) :: ddobsini    ! Obs. ini time in YYYYMMDD.HHMMSS 
    86       REAL(KIND=dp), INTENT(IN) :: ddobsend    ! Obs. end time in YYYYMMDD.HHMMSS 
     74      INTEGER, INTENT(IN) :: kextr      ! Number of extra fields for each var 
     75      INTEGER, INTENT(IN) :: kstp       ! Ocean time-step index 
     76      LOGICAL, INTENT(IN) :: ldvar1     ! Observed variables switches 
     77      LOGICAL, INTENT(IN) :: ldvar2 
     78      LOGICAL, INTENT(IN) :: ldignmis   ! Ignore missing files 
     79      LOGICAL, INTENT(IN) :: ldsatt     ! Compute salinity at all temperature points 
     80      LOGICAL, INTENT(IN) :: ldmod      ! Initialize model from input data 
     81      REAL(dp), INTENT(IN) :: ddobsini  ! Obs. ini time in YYYYMMDD.HHMMSS 
     82      REAL(dp), INTENT(IN) :: ddobsend  ! Obs. end time in YYYYMMDD.HHMMSS 
    8783      INTEGER, DIMENSION(imaxavtypes), OPTIONAL :: & 
    88          & kdailyavtypes 
     84         & kdailyavtypes                ! Types of daily average observations 
    8985 
    9086      !! * Local declarations 
    91       CHARACTER(LEN=15), PARAMETER :: cpname='obs_rea_pro_dri' 
     87      CHARACTER(LEN=15), PARAMETER :: cpname='obs_rea_prof' 
     88      CHARACTER(len=8) :: clrefdate 
     89      CHARACTER(len=6), DIMENSION(:), ALLOCATABLE :: clvars 
    9290      INTEGER :: jvar 
    9391      INTEGER :: ji 
     
    105103      INTEGER :: imin 
    106104      INTEGER :: isec 
     105      INTEGER :: iprof 
     106      INTEGER :: iproftot 
     107      INTEGER :: ivar1t0 
     108      INTEGER :: ivar2t0 
     109      INTEGER :: ivar1t 
     110      INTEGER :: ivar2t 
     111      INTEGER :: ip3dt 
     112      INTEGER :: ios 
     113      INTEGER :: ioserrcount 
     114      INTEGER :: ivar1tmpp 
     115      INTEGER :: ivar2tmpp 
     116      INTEGER :: ip3dtmpp 
     117      INTEGER :: itype 
    107118      INTEGER, DIMENSION(knumfiles) :: & 
    108119         & irefdate 
    109120      INTEGER, DIMENSION(ntyp1770+1) :: & 
    110          & itypt,    & 
    111          & ityptmpp, & 
    112          & ityps,    & 
    113          & itypsmpp  
    114       INTEGER :: it3dtmpp 
    115       INTEGER :: is3dtmpp 
    116       INTEGER :: ip3dtmpp 
     121         & itypvar1,    & 
     122         & itypvar1mpp, & 
     123         & itypvar2,    & 
     124         & itypvar2mpp  
    117125      INTEGER, DIMENSION(:), ALLOCATABLE :: & 
    118          & iobsi,    & 
    119          & iobsj,    & 
    120          & iproc,    & 
     126         & iobsi1,    & 
     127         & iobsj1,    & 
     128         & iproc1,    & 
     129         & iobsi2,    & 
     130         & iobsj2,    & 
     131         & iproc2,    & 
    121132         & iindx,    & 
    122133         & ifileidx, & 
    123134         & iprofidx 
    124       INTEGER :: itype 
    125135      INTEGER, DIMENSION(imaxavtypes) :: & 
    126136         & idailyavtypes 
     137      INTEGER, DIMENSION(kvars) :: & 
     138         & iv3dt 
    127139      REAL(wp), DIMENSION(:), ALLOCATABLE :: & 
    128140         & zphi, & 
    129141         & zlam 
    130       real(wp), DIMENSION(:), ALLOCATABLE :: & 
     142      REAL(wp), DIMENSION(:), ALLOCATABLE :: & 
    131143         & zdat 
     144      REAL(wp), DIMENSION(knumfiles) :: & 
     145         & djulini, & 
     146         & djulend 
    132147      LOGICAL :: llvalprof 
     148      LOGICAL :: lldavtimset 
    133149      TYPE(obfbdata), POINTER, DIMENSION(:) :: & 
    134150         & inpfiles 
    135       real(wp), DIMENSION(knumfiles) :: & 
    136          & djulini, & 
    137          & djulend 
    138       INTEGER :: iprof 
    139       INTEGER :: iproftot 
    140       INTEGER :: it3dt0 
    141       INTEGER :: is3dt0 
    142       INTEGER :: it3dt 
    143       INTEGER :: is3dt 
    144       INTEGER :: ip3dt 
    145       INTEGER :: ios 
    146       INTEGER :: ioserrcount 
    147       INTEGER, DIMENSION(kvars) :: & 
    148          & iv3dt 
    149       CHARACTER(len=8) :: cl_refdate 
    150     
     151 
    151152      ! Local initialization 
    152153      iprof = 0 
    153       it3dt0 = 0 
    154       is3dt0 = 0 
     154      ivar1t0 = 0 
     155      ivar2t0 = 0 
    155156      ip3dt = 0 
    156157 
    157158      ! Daily average types 
     159      lldavtimset = .FALSE. 
    158160      IF ( PRESENT(kdailyavtypes) ) THEN 
    159161         idailyavtypes(:) = kdailyavtypes(:) 
     162         IF ( ANY (idailyavtypes(:) /= -1) ) lldavtimset = .TRUE. 
    160163      ELSE 
    161164         idailyavtypes(:) = -1 
     
    163166 
    164167      !----------------------------------------------------------------------- 
    165       ! Check data the model part is just with feedback data files 
    166       !----------------------------------------------------------------------- 
    167       IF ( ldmod .AND. ( kformat /= 0 ) ) THEN 
    168          CALL ctl_stop( 'Model can only be read from feedback data' ) 
    169          RETURN 
    170       ENDIF 
    171  
    172       !----------------------------------------------------------------------- 
    173168      ! Count the number of files needed and allocate the obfbdata type 
    174169      !----------------------------------------------------------------------- 
    175        
     170 
    176171      inobf = knumfiles 
    177        
     172 
    178173      ALLOCATE( inpfiles(inobf) ) 
    179174 
    180175      prof_files : DO jj = 1, inobf 
    181            
     176 
    182177         !--------------------------------------------------------------------- 
    183178         ! Prints 
     
    186181            WRITE(numout,*) 
    187182            WRITE(numout,*) ' obs_rea_pro_dri : Reading from file = ', & 
    188                & TRIM( TRIM( cfilenames(jj) ) ) 
     183               & TRIM( TRIM( cdfilenames(jj) ) ) 
    189184            WRITE(numout,*) ' ~~~~~~~~~~~~~~~' 
    190185            WRITE(numout,*) 
     
    194189         !  Initialization: Open file and get dimensions only 
    195190         !--------------------------------------------------------------------- 
    196           
    197          iflag = nf90_open( TRIM( TRIM( cfilenames(jj) ) ), nf90_nowrite, & 
     191 
     192         iflag = nf90_open( TRIM( cdfilenames(jj) ), nf90_nowrite, & 
    198193            &                      i_file_id ) 
    199           
     194 
    200195         IF ( iflag /= nf90_noerr ) THEN 
    201196 
    202197            IF ( ldignmis ) THEN 
    203198               inpfiles(jj)%nobs = 0 
    204                CALL ctl_warn( 'File ' // TRIM( TRIM( cfilenames(jj) ) ) // & 
     199               CALL ctl_warn( 'File ' // TRIM( cdfilenames(jj) ) // & 
    205200                  &           ' not found' ) 
    206201            ELSE  
    207                CALL ctl_stop( 'File ' // TRIM( TRIM( cfilenames(jj) ) ) // & 
     202               CALL ctl_stop( 'File ' // TRIM( cdfilenames(jj) ) // & 
    208203                  &           ' not found' ) 
    209204            ENDIF 
    210205 
    211206         ELSE  
    212              
     207 
    213208            !------------------------------------------------------------------ 
    214             !  Close the file since it is opened in read_proffile 
     209            !  Close the file since it is opened in read_obfbdata 
    215210            !------------------------------------------------------------------ 
    216              
     211 
    217212            iflag = nf90_close( i_file_id ) 
    218213 
     
    220215            !  Read the profile file into inpfiles 
    221216            !------------------------------------------------------------------ 
    222             IF ( kformat == 0 ) THEN 
    223                CALL init_obfbdata( inpfiles(jj) ) 
    224                IF(lwp) THEN 
    225                   WRITE(numout,*) 
    226                   WRITE(numout,*)'Reading from feedback file :', & 
    227                      &           TRIM( cfilenames(jj) ) 
    228                ENDIF 
    229                CALL read_obfbdata( TRIM( cfilenames(jj) ), inpfiles(jj), & 
    230                   &                ldgrid = .TRUE. ) 
    231                IF ( inpfiles(jj)%nvar < 2 ) THEN 
    232                   CALL ctl_stop( 'Feedback format error' ) 
    233                   RETURN 
    234                ENDIF 
    235                IF ( TRIM(inpfiles(jj)%cname(1)) /= 'POTM' ) THEN 
    236                   CALL ctl_stop( 'Feedback format error' ) 
    237                   RETURN 
    238                ENDIF 
    239                IF ( TRIM(inpfiles(jj)%cname(2)) /= 'PSAL' ) THEN 
    240                   CALL ctl_stop( 'Feedback format error' ) 
    241                   RETURN 
    242                ENDIF 
    243                IF ( ldmod .AND. ( inpfiles(jj)%nadd == 0 ) ) THEN 
    244                   CALL ctl_stop( 'Model not in input data' ) 
    245                   RETURN 
    246                ENDIF 
    247             ELSEIF ( kformat == 1 ) THEN 
    248                CALL read_enactfile( TRIM( cfilenames(jj) ), inpfiles(jj), & 
    249                   &                 numout, lwp, .TRUE. ) 
    250             ELSEIF ( kformat == 2 ) THEN 
    251                CALL read_coriofile( TRIM( cfilenames(jj) ), inpfiles(jj), & 
    252                   &                 numout, lwp, .TRUE. ) 
     217            CALL init_obfbdata( inpfiles(jj) ) 
     218            CALL read_obfbdata( TRIM( cdfilenames(jj) ), inpfiles(jj), & 
     219               &                ldgrid = .TRUE. ) 
     220 
     221            IF ( inpfiles(jj)%nvar < 2 ) THEN 
     222               CALL ctl_stop( 'Feedback format error: ', & 
     223                  &           ' less than 2 vars in profile file' ) 
     224            ENDIF 
     225 
     226            IF ( ldmod .AND. ( inpfiles(jj)%nadd == 0 ) ) THEN 
     227               CALL ctl_stop( 'Model not in input data' ) 
     228            ENDIF 
     229 
     230            IF ( jj == 1 ) THEN 
     231               ALLOCATE( clvars( inpfiles(jj)%nvar ) ) 
     232               DO ji = 1, inpfiles(jj)%nvar 
     233                 clvars(ji) = inpfiles(jj)%cname(ji) 
     234               END DO 
    253235            ELSE 
    254                CALL ctl_stop( 'File format unknown' ) 
    255             ENDIF 
    256              
     236               DO ji = 1, inpfiles(jj)%nvar 
     237                  IF ( inpfiles(jj)%cname(ji) /= clvars(ji) ) THEN 
     238                     CALL ctl_stop( 'Feedback file variables not consistent', & 
     239                        &           ' with previous files for this type' ) 
     240                  ENDIF 
     241               END DO 
     242            ENDIF 
     243 
    257244            !------------------------------------------------------------------ 
    258245            !  Change longitude (-180,180) 
     
    272259            !  Calculate the date  (change eventually) 
    273260            !------------------------------------------------------------------ 
    274             cl_refdate=inpfiles(jj)%cdjuldref(1:8) 
    275             READ(cl_refdate,'(I8)') irefdate(jj) 
    276              
     261            clrefdate=inpfiles(jj)%cdjuldref(1:8) 
     262            READ(clrefdate,'(I8)') irefdate(jj) 
     263 
    277264            CALL ddatetoymdhms( ddobsini, iyea, imon, iday, ihou, imin, isec ) 
    278265            CALL greg2jul( isec, imin, ihou, iday, imon, iyea, djulini(jj), & 
     
    283270 
    284271            ioserrcount=0 
    285             IF ( ldavtimset ) THEN 
     272            IF ( lldavtimset ) THEN 
     273 
     274               IF ( ANY ( idailyavtypes(:) /= -1 ) .AND. lwp) THEN 
     275                  WRITE(numout,*)' Resetting time of daily averaged', & 
     276                     &           ' observations to the end of the day' 
     277               ENDIF 
     278 
    286279               DO ji = 1, inpfiles(jj)%nobs 
    287                   !  
    288                   !  for daily averaged data for example 
    289                   !  MRB data (itype==820) force the time 
    290                   !  to be the  end of the day 
    291                   ! 
    292280                  READ( inpfiles(jj)%cdtyp(ji), '(I4)', IOSTAT = ios, ERR = 900 ) itype 
    293281900               IF ( ios /= 0 ) THEN 
    294                      itype = 0         ! Set type to zero if there is a problem in the string conversion 
    295                   ENDIF 
    296                   IF ( ANY (idailyavtypes == itype ) ) THEN 
    297                      inpfiles(jj)%ptim(ji) = & 
    298                      & INT(inpfiles(jj)%ptim(ji)) + 1 
    299                   ENDIF 
     282                     ! Set type to zero if there is a problem in the string conversion 
     283                     itype = 0 
     284                  ENDIF 
     285 
     286                  IF ( ANY ( idailyavtypes(:) == itype ) ) THEN 
     287                  !  for daily averaged data force the time 
     288                  !  to be the last time-step of the day, but still within the day. 
     289                     IF ( inpfiles(jj)%ptim(ji) >= 0. ) THEN 
     290                        inpfiles(jj)%ptim(ji) = & 
     291                           & INT(inpfiles(jj)%ptim(ji)) + 0.9999 
     292                     ELSE 
     293                        inpfiles(jj)%ptim(ji) = & 
     294                           & INT(inpfiles(jj)%ptim(ji)) - 0.0001 
     295                     ENDIF 
     296                  ENDIF 
     297 
    300298               END DO 
    301             ENDIF 
    302              
     299 
     300            ENDIF 
     301 
    303302            IF ( inpfiles(jj)%nobs > 0 ) THEN 
    304                inpfiles(jj)%iproc = -1 
    305                inpfiles(jj)%iobsi = -1 
    306                inpfiles(jj)%iobsj = -1 
     303               inpfiles(jj)%iproc(:,:) = -1 
     304               inpfiles(jj)%iobsi(:,:) = -1 
     305               inpfiles(jj)%iobsj(:,:) = -1 
    307306            ENDIF 
    308307            inowin = 0 
     
    318317            ALLOCATE( zlam(inowin)  ) 
    319318            ALLOCATE( zphi(inowin)  ) 
    320             ALLOCATE( iobsi(inowin) ) 
    321             ALLOCATE( iobsj(inowin) ) 
    322             ALLOCATE( iproc(inowin) ) 
     319            ALLOCATE( iobsi1(inowin) ) 
     320            ALLOCATE( iobsj1(inowin) ) 
     321            ALLOCATE( iproc1(inowin) ) 
     322            ALLOCATE( iobsi2(inowin) ) 
     323            ALLOCATE( iobsj2(inowin) ) 
     324            ALLOCATE( iproc2(inowin) ) 
    323325            inowin = 0 
    324326            DO ji = 1, inpfiles(jj)%nobs 
     
    334336            END DO 
    335337 
    336             CALL obs_grid_search( inowin, zlam, zphi, iobsi, iobsj, iproc, 'T' ) 
     338            IF ( TRIM( inpfiles(jj)%cname(1) ) == 'POTM' ) THEN 
     339               CALL obs_grid_search( inowin, zlam, zphi, iobsi1, iobsj1, & 
     340                  &                  iproc1, 'T' ) 
     341               iobsi2(:) = iobsi1(:) 
     342               iobsj2(:) = iobsj1(:) 
     343               iproc2(:) = iproc1(:) 
     344            ELSEIF ( TRIM( inpfiles(jj)%cname(1) ) == 'UVEL' ) THEN 
     345               CALL obs_grid_search( inowin, zlam, zphi, iobsi1, iobsj1, & 
     346                  &                  iproc1, 'U' ) 
     347               CALL obs_grid_search( inowin, zlam, zphi, iobsi2, iobsj2, & 
     348                  &                  iproc2, 'V' ) 
     349            ENDIF 
    337350 
    338351            inowin = 0 
     
    344357                  & ( inpfiles(jj)%ptim(ji) <= djulend(jj) )       ) THEN 
    345358                  inowin = inowin + 1 
    346                   inpfiles(jj)%iproc(ji,1) = iproc(inowin) 
    347                   inpfiles(jj)%iobsi(ji,1) = iobsi(inowin) 
    348                   inpfiles(jj)%iobsj(ji,1) = iobsj(inowin) 
     359                  inpfiles(jj)%iproc(ji,1) = iproc1(inowin) 
     360                  inpfiles(jj)%iobsi(ji,1) = iobsi1(inowin) 
     361                  inpfiles(jj)%iobsj(ji,1) = iobsj1(inowin) 
     362                  inpfiles(jj)%iproc(ji,2) = iproc2(inowin) 
     363                  inpfiles(jj)%iobsi(ji,2) = iobsi2(inowin) 
     364                  inpfiles(jj)%iobsj(ji,2) = iobsj2(inowin) 
     365                  IF ( inpfiles(jj)%iproc(ji,1) /= & 
     366                     & inpfiles(jj)%iproc(ji,2) ) THEN 
     367                     CALL ctl_stop( 'Error in obs_read_prof:', & 
     368                        & 'var1 and var2 observation on different processors') 
     369                  ENDIF 
    349370               ENDIF 
    350371            END DO 
    351             DEALLOCATE( zlam, zphi, iobsi, iobsj, iproc ) 
     372            DEALLOCATE( zlam, zphi, iobsi1, iobsj1, iproc1, iobsi2, iobsj2, iproc2 ) 
    352373 
    353374            DO ji = 1, inpfiles(jj)%nobs 
     
    363384                  ENDIF 
    364385                  llvalprof = .FALSE. 
    365                   IF ( ldt3d ) THEN 
     386                  IF ( ldvar1 ) THEN 
    366387                     loop_t_count : DO ij = 1,inpfiles(jj)%nlev 
    367388                        IF ( inpfiles(jj)%pdep(ij,ji) >= 6000. ) & 
     
    369390                        IF ( ( inpfiles(jj)%ivlqc(ij,ji,1) <= 2 ) .AND. & 
    370391                           & ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) THEN 
    371                            it3dt0 = it3dt0 + 1 
     392                           ivar1t0 = ivar1t0 + 1 
    372393                        ENDIF 
    373394                     END DO loop_t_count 
    374395                  ENDIF 
    375                   IF ( lds3d ) THEN 
     396                  IF ( ldvar2 ) THEN 
    376397                     loop_s_count : DO ij = 1,inpfiles(jj)%nlev 
    377398                        IF ( inpfiles(jj)%pdep(ij,ji) >= 6000. ) & 
     
    379400                        IF ( ( inpfiles(jj)%ivlqc(ij,ji,2) <= 2 ) .AND. & 
    380401                           & ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) THEN 
    381                            is3dt0 = is3dt0 + 1 
     402                           ivar2t0 = ivar2t0 + 1 
    382403                        ENDIF 
    383404                     END DO loop_s_count 
     
    388409                     IF ( ( ( ( inpfiles(jj)%ivlqc(ij,ji,1) <= 2 ) .AND. & 
    389410                        &     ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) .AND. & 
    390                         &     ldt3d ) .OR. & 
     411                        &     ldvar1 ) .OR. & 
    391412                        & ( ( ( inpfiles(jj)%ivlqc(ij,ji,2) <= 2 ) .AND. & 
    392413                        &     ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) .AND. & 
    393                         &     lds3d ) ) THEN 
     414                        &     ldvar2 ) ) THEN 
    394415                        ip3dt = ip3dt + 1 
    395416                        llvalprof = .TRUE. 
     
    405426 
    406427      END DO prof_files 
    407        
     428 
    408429      !----------------------------------------------------------------------- 
    409430      ! Get the time ordered indices of the input data 
     
    446467         &               zdat,     & 
    447468         &               iindx   ) 
    448        
     469 
    449470      iv3dt(:) = -1 
    450471      IF (ldsatt) THEN 
     
    452473         iv3dt(2) = ip3dt 
    453474      ELSE 
    454          iv3dt(1) = it3dt0 
    455          iv3dt(2) = is3dt0 
     475         iv3dt(1) = ivar1t0 
     476         iv3dt(2) = ivar2t0 
    456477      ENDIF 
    457478      CALL obs_prof_alloc( profdata, kvars, kextr, iprof, iv3dt, & 
    458479         &                 kstp, jpi, jpj, jpk ) 
    459        
     480 
    460481      ! * Read obs/positions, QC, all variable and assign to profdata 
    461482 
    462483      profdata%nprof     = 0 
    463484      profdata%nvprot(:) = 0 
    464  
     485      profdata%cvars(:)  = clvars(:) 
    465486      iprof = 0 
    466487 
    467488      ip3dt = 0 
    468       it3dt = 0 
    469       is3dt = 0 
    470       itypt   (:) = 0 
    471       ityptmpp(:) = 0 
    472        
    473       ityps   (:) = 0 
    474       itypsmpp(:) = 0 
    475        
    476       ioserrcount = 0       
     489      ivar1t = 0 
     490      ivar2t = 0 
     491      itypvar1   (:) = 0 
     492      itypvar1mpp(:) = 0 
     493 
     494      itypvar2   (:) = 0 
     495      itypvar2mpp(:) = 0 
     496 
     497      ioserrcount = 0 
    477498      DO jk = 1, iproftot 
    478           
     499 
    479500         jj = ifileidx(iindx(jk)) 
    480501         ji = iprofidx(iindx(jk)) 
     
    486507         IF ( ( inpfiles(jj)%ptim(ji) >  djulini(jj) ) .AND.  & 
    487508            & ( inpfiles(jj)%ptim(ji) <= djulend(jj) ) ) THEN 
    488              
     509 
    489510            IF ( nproc == 0 ) THEN 
    490511               IF ( inpfiles(jj)%iproc(ji,1) >  nproc ) CYCLE 
     
    492513               IF ( inpfiles(jj)%iproc(ji,1) /= nproc ) CYCLE 
    493514            ENDIF 
    494              
     515 
    495516            llvalprof = .FALSE. 
    496517 
     
    501522 
    502523            loop_prof : DO ij = 1, inpfiles(jj)%nlev 
    503                 
     524 
    504525               IF ( inpfiles(jj)%pdep(ij,ji) >= 6000. ) & 
    505526                  & CYCLE 
    506                 
     527 
    507528               IF ( ( inpfiles(jj)%ivlqc(ij,ji,1) <= 2 ) .AND. & 
    508529                  & ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) THEN 
    509                    
     530 
    510531                  llvalprof = .TRUE.  
    511532                  EXIT loop_prof 
    512                    
     533 
    513534               ENDIF 
    514                 
     535 
    515536               IF ( ( inpfiles(jj)%ivlqc(ij,ji,2) <= 2 ) .AND. & 
    516537                  & ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) THEN 
    517                    
     538 
    518539                  llvalprof = .TRUE.  
    519540                  EXIT loop_prof 
    520                    
     541 
    521542               ENDIF 
    522                 
     543 
    523544            END DO loop_prof 
    524              
     545 
    525546            ! Set profile information 
    526              
     547 
    527548            IF ( llvalprof ) THEN 
    528                 
     549 
    529550               iprof = iprof + 1 
    530551 
     
    545566               profdata%nhou(iprof) = ihou 
    546567               profdata%nmin(iprof) = imin 
    547                 
     568 
    548569               ! Profile space coordinates 
    549570               profdata%rlam(iprof) = inpfiles(jj)%plam(ji) 
     
    551572 
    552573               ! Coordinate search parameters 
    553                profdata%mi  (iprof,:) = inpfiles(jj)%iobsi(ji,1) 
    554                profdata%mj  (iprof,:) = inpfiles(jj)%iobsj(ji,1) 
    555                 
     574               profdata%mi  (iprof,1) = inpfiles(jj)%iobsi(ji,1) 
     575               profdata%mj  (iprof,1) = inpfiles(jj)%iobsj(ji,1) 
     576               profdata%mi  (iprof,2) = inpfiles(jj)%iobsi(ji,2) 
     577               profdata%mj  (iprof,2) = inpfiles(jj)%iobsj(ji,2) 
     578 
    556579               ! Profile WMO number 
    557580               profdata%cwmo(iprof) = inpfiles(jj)%cdwmo(ji) 
    558                 
     581 
    559582               ! Instrument type 
    560583               READ( inpfiles(jj)%cdtyp(ji), '(I4)', IOSTAT = ios, ERR = 901 ) itype 
     
    564587                  itype = 0 
    565588               ENDIF 
    566                 
     589 
    567590               profdata%ntyp(iprof) = itype 
    568                 
     591 
    569592               ! QC stuff 
    570593 
     
    585608               profdata%nqc(iprof)  = 0 !TODO 
    586609 
    587                loop_p : DO ij = 1, inpfiles(jj)%nlev             
    588                    
     610               loop_p : DO ij = 1, inpfiles(jj)%nlev 
     611 
    589612                  IF ( inpfiles(jj)%pdep(ij,ji) >= 6000. ) & 
    590613                     & CYCLE 
     
    594617                     IF ( ( ( ( inpfiles(jj)%ivlqc(ij,ji,1) <= 2 ) .AND. & 
    595618                        &     ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) .AND. & 
    596                         &     ldt3d ) .OR. & 
     619                        &     ldvar1 ) .OR. & 
    597620                        & ( ( ( inpfiles(jj)%ivlqc(ij,ji,2) <= 2 ) .AND. & 
    598621                        &     ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) .AND. & 
    599                         &     lds3d ) ) THEN 
     622                        &     ldvar2 ) ) THEN 
    600623                        ip3dt = ip3dt + 1 
    601624                     ELSE 
    602625                        CYCLE 
    603626                     ENDIF 
    604                       
     627 
    605628                  ENDIF 
    606629 
    607630                  IF ( ( ( ( inpfiles(jj)%ivlqc(ij,ji,1) <= 2 ) .AND. & 
    608631                     &     ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) .AND. & 
    609                      &       ldt3d ) .OR. ldsatt ) THEN 
    610                       
     632                     &       ldvar1 ) .OR. ldsatt ) THEN 
     633 
    611634                     IF (ldsatt) THEN 
    612635 
    613                         it3dt = ip3dt 
     636                        ivar1t = ip3dt 
    614637 
    615638                     ELSE 
    616639 
    617                         it3dt = it3dt + 1 
    618                          
     640                        ivar1t = ivar1t + 1 
     641 
    619642                     ENDIF 
    620643 
    621                      ! Depth of T observation 
    622                      profdata%var(1)%vdep(it3dt) = & 
     644                     ! Depth of var1 observation 
     645                     profdata%var(1)%vdep(ivar1t) = & 
    623646                        &                inpfiles(jj)%pdep(ij,ji) 
    624                       
    625                      ! Depth of T observation QC 
    626                      profdata%var(1)%idqc(it3dt) = & 
     647 
     648                     ! Depth of var1 observation QC 
     649                     profdata%var(1)%idqc(ivar1t) = & 
    627650                        &                inpfiles(jj)%idqc(ij,ji) 
    628                       
    629                      ! Depth of T observation QC flags 
    630                      profdata%var(1)%idqcf(:,it3dt) = & 
     651 
     652                     ! Depth of var1 observation QC flags 
     653                     profdata%var(1)%idqcf(:,ivar1t) = & 
    631654                        &                inpfiles(jj)%idqcf(:,ij,ji) 
    632                       
     655 
    633656                     ! Profile index 
    634                      profdata%var(1)%nvpidx(it3dt) = iprof 
    635                       
     657                     profdata%var(1)%nvpidx(ivar1t) = iprof 
     658 
    636659                     ! Vertical index in original profile 
    637                      profdata%var(1)%nvlidx(it3dt) = ij 
    638  
    639                      ! Profile potential T value 
     660                     profdata%var(1)%nvlidx(ivar1t) = ij 
     661 
     662                     ! Profile var1 value 
    640663                     IF ( ( inpfiles(jj)%ivlqc(ij,ji,1) <= 2 ) .AND. & 
    641664                        & ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) THEN 
    642                         profdata%var(1)%vobs(it3dt) = & 
     665                        profdata%var(1)%vobs(ivar1t) = & 
    643666                           &                inpfiles(jj)%pob(ij,ji,1) 
    644667                        IF ( ldmod ) THEN 
    645                            profdata%var(1)%vmod(it3dt) = & 
     668                           profdata%var(1)%vmod(ivar1t) = & 
    646669                              &                inpfiles(jj)%padd(ij,ji,1,1) 
    647670                        ENDIF 
    648                         ! Count number of profile T data as function of type 
    649                         itypt( profdata%ntyp(iprof) + 1 ) = & 
    650                            & itypt( profdata%ntyp(iprof) + 1 ) + 1 
     671                        ! Count number of profile var1 data as function of type 
     672                        itypvar1( profdata%ntyp(iprof) + 1 ) = & 
     673                           & itypvar1( profdata%ntyp(iprof) + 1 ) + 1 
    651674                     ELSE 
    652                         profdata%var(1)%vobs(it3dt) = fbrmdi 
     675                        profdata%var(1)%vobs(ivar1t) = fbrmdi 
    653676                     ENDIF 
    654677 
    655                      ! Profile T qc 
    656                      profdata%var(1)%nvqc(it3dt) = & 
     678                     ! Profile var1 qc 
     679                     profdata%var(1)%nvqc(ivar1t) = & 
    657680                        & inpfiles(jj)%ivlqc(ij,ji,1) 
    658681 
    659                      ! Profile T qc flags 
    660                      profdata%var(1)%nvqcf(:,it3dt) = & 
     682                     ! Profile var1 qc flags 
     683                     profdata%var(1)%nvqcf(:,ivar1t) = & 
    661684                        & inpfiles(jj)%ivlqcf(:,ij,ji,1) 
    662685 
    663686                     ! Profile insitu T value 
    664                      profdata%var(1)%vext(it3dt,1) = & 
    665                         &                inpfiles(jj)%pext(ij,ji,1) 
    666                       
    667                   ENDIF 
    668                    
     687                     IF ( TRIM( inpfiles(jj)%cname(1) ) == 'POTM' ) THEN 
     688                        profdata%var(1)%vext(ivar1t,1) = & 
     689                           &                inpfiles(jj)%pext(ij,ji,1) 
     690                     ENDIF 
     691 
     692                  ENDIF 
     693 
    669694                  IF ( ( ( ( inpfiles(jj)%ivlqc(ij,ji,2) <= 2 ) .AND. & 
    670695                     &   ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) .AND. & 
    671                      &   lds3d ) .OR. ldsatt ) THEN 
    672                       
     696                     &   ldvar2 ) .OR. ldsatt ) THEN 
     697 
    673698                     IF (ldsatt) THEN 
    674699 
    675                         is3dt = ip3dt 
     700                        ivar2t = ip3dt 
    676701 
    677702                     ELSE 
    678703 
    679                         is3dt = is3dt + 1 
    680                          
     704                        ivar2t = ivar2t + 1 
     705 
    681706                     ENDIF 
    682707 
    683                      ! Depth of S observation 
    684                      profdata%var(2)%vdep(is3dt) = & 
     708                     ! Depth of var2 observation 
     709                     profdata%var(2)%vdep(ivar2t) = & 
    685710                        &                inpfiles(jj)%pdep(ij,ji) 
    686                       
    687                      ! Depth of S observation QC 
    688                      profdata%var(2)%idqc(is3dt) = & 
     711 
     712                     ! Depth of var2 observation QC 
     713                     profdata%var(2)%idqc(ivar2t) = & 
    689714                        &                inpfiles(jj)%idqc(ij,ji) 
    690                       
    691                      ! Depth of S observation QC flags 
    692                      profdata%var(2)%idqcf(:,is3dt) = & 
     715 
     716                     ! Depth of var2 observation QC flags 
     717                     profdata%var(2)%idqcf(:,ivar2t) = & 
    693718                        &                inpfiles(jj)%idqcf(:,ij,ji) 
    694                       
     719 
    695720                     ! Profile index 
    696                      profdata%var(2)%nvpidx(is3dt) = iprof 
    697                       
     721                     profdata%var(2)%nvpidx(ivar2t) = iprof 
     722 
    698723                     ! Vertical index in original profile 
    699                      profdata%var(2)%nvlidx(is3dt) = ij 
    700  
    701                      ! Profile S value 
     724                     profdata%var(2)%nvlidx(ivar2t) = ij 
     725 
     726                     ! Profile var2 value 
    702727                     IF ( ( inpfiles(jj)%ivlqc(ij,ji,2) <= 2 ) .AND. & 
    703728                        & ( inpfiles(jj)%idqc(ij,ji) <= 2 ) ) THEN 
    704                         profdata%var(2)%vobs(is3dt) = & 
     729                        profdata%var(2)%vobs(ivar2t) = & 
    705730                           &                inpfiles(jj)%pob(ij,ji,2) 
    706731                        IF ( ldmod ) THEN 
    707                            profdata%var(2)%vmod(is3dt) = & 
     732                           profdata%var(2)%vmod(ivar2t) = & 
    708733                              &                inpfiles(jj)%padd(ij,ji,1,2) 
    709734                        ENDIF 
    710                         ! Count number of profile S data as function of type 
    711                         ityps( profdata%ntyp(iprof) + 1 ) = & 
    712                            & ityps( profdata%ntyp(iprof) + 1 ) + 1 
     735                        ! Count number of profile var2 data as function of type 
     736                        itypvar2( profdata%ntyp(iprof) + 1 ) = & 
     737                           & itypvar2( profdata%ntyp(iprof) + 1 ) + 1 
    713738                     ELSE 
    714                         profdata%var(2)%vobs(is3dt) = fbrmdi 
     739                        profdata%var(2)%vobs(ivar2t) = fbrmdi 
    715740                     ENDIF 
    716                       
    717                      ! Profile S qc 
    718                      profdata%var(2)%nvqc(is3dt) = & 
     741 
     742                     ! Profile var2 qc 
     743                     profdata%var(2)%nvqc(ivar2t) = & 
    719744                        & inpfiles(jj)%ivlqc(ij,ji,2) 
    720745 
    721                      ! Profile S qc flags 
    722                      profdata%var(2)%nvqcf(:,is3dt) = & 
     746                     ! Profile var2 qc flags 
     747                     profdata%var(2)%nvqcf(:,ivar2t) = & 
    723748                        & inpfiles(jj)%ivlqcf(:,ij,ji,2) 
    724749 
    725750                  ENDIF 
    726              
     751 
    727752               END DO loop_p 
    728753 
     
    736761      ! Sum up over processors 
    737762      !----------------------------------------------------------------------- 
    738        
    739       CALL obs_mpp_sum_integer ( it3dt0, it3dtmpp ) 
    740       CALL obs_mpp_sum_integer ( is3dt0, is3dtmpp ) 
    741       CALL obs_mpp_sum_integer ( ip3dt, ip3dtmpp ) 
    742        
    743       CALL obs_mpp_sum_integers( itypt, ityptmpp, ntyp1770 + 1 ) 
    744       CALL obs_mpp_sum_integers( ityps, itypsmpp, ntyp1770 + 1 ) 
    745        
     763 
     764      CALL obs_mpp_sum_integer ( ivar1t0, ivar1tmpp ) 
     765      CALL obs_mpp_sum_integer ( ivar2t0, ivar2tmpp ) 
     766      CALL obs_mpp_sum_integer ( ip3dt,   ip3dtmpp ) 
     767 
     768      CALL obs_mpp_sum_integers( itypvar1, itypvar1mpp, ntyp1770 + 1 ) 
     769      CALL obs_mpp_sum_integers( itypvar2, itypvar2mpp, ntyp1770 + 1 ) 
     770 
    746771      !----------------------------------------------------------------------- 
    747772      ! Output number of observations. 
     
    749774      IF(lwp) THEN 
    750775         WRITE(numout,*)  
    751          WRITE(numout,'(1X,A)') 'Profile data' 
     776         WRITE(numout,'(A)') ' Profile data' 
    752777         WRITE(numout,'(1X,A)') '------------' 
    753778         WRITE(numout,*)  
    754          WRITE(numout,'(1X,A)') 'Profile T data' 
    755          WRITE(numout,'(1X,A)') '--------------' 
     779         WRITE(numout,'(1X,A)') 'Profile data, '//TRIM( profdata%cvars(1) ) 
     780         WRITE(numout,'(1X,A)') '------------------------' 
    756781         DO ji = 0, ntyp1770 
    757             IF ( ityptmpp(ji+1) > 0 ) THEN 
     782            IF ( itypvar1mpp(ji+1) > 0 ) THEN 
    758783               WRITE(numout,'(1X,A3,1X,A48,A3,I8)') ctypshort(ji), & 
    759784                  & cwmonam1770(ji)(1:52),' = ', & 
    760                   & ityptmpp(ji+1) 
     785                  & itypvar1mpp(ji+1) 
    761786            ENDIF 
    762787         END DO 
     
    764789            & '---------------------------------------------------------------' 
    765790         WRITE(numout,'(1X,A55,I8)') & 
    766             & 'Total profile T data                                 = ',& 
    767             & it3dtmpp 
     791            & 'Total profile data for variable '//TRIM( profdata%cvars(1) )// & 
     792            & '             = ', ivar1tmpp 
    768793         WRITE(numout,'(1X,A)') & 
    769794            & '---------------------------------------------------------------' 
    770795         WRITE(numout,*)  
    771          WRITE(numout,'(1X,A)') 'Profile S data' 
    772          WRITE(numout,'(1X,A)') '--------------' 
     796         WRITE(numout,'(1X,A)') 'Profile data, '//TRIM( profdata%cvars(2) ) 
     797         WRITE(numout,'(1X,A)') '------------------------' 
    773798         DO ji = 0, ntyp1770 
    774             IF ( itypsmpp(ji+1) > 0 ) THEN 
     799            IF ( itypvar2mpp(ji+1) > 0 ) THEN 
    775800               WRITE(numout,'(1X,A3,1X,A48,A3,I8)') ctypshort(ji), & 
    776801                  & cwmonam1770(ji)(1:52),' = ', & 
    777                   & itypsmpp(ji+1) 
     802                  & itypvar2mpp(ji+1) 
    778803            ENDIF 
    779804         END DO 
     
    781806            & '---------------------------------------------------------------' 
    782807         WRITE(numout,'(1X,A55,I8)') & 
    783             & 'Total profile S data                                 = ',& 
    784             & is3dtmpp 
     808            & 'Total profile data for variable '//TRIM( profdata%cvars(2) )// & 
     809            & '             = ', ivar2tmpp 
    785810         WRITE(numout,'(1X,A)') & 
    786811            & '---------------------------------------------------------------' 
    787812         WRITE(numout,*)  
    788813      ENDIF 
    789        
     814 
    790815      IF (ldsatt) THEN 
    791816         profdata%nvprot(1)    = ip3dt 
     
    794819         profdata%nvprotmpp(2) = ip3dtmpp 
    795820      ELSE 
    796          profdata%nvprot(1)    = it3dt 
    797          profdata%nvprot(2)    = is3dt 
    798          profdata%nvprotmpp(1) = it3dtmpp 
    799          profdata%nvprotmpp(2) = is3dtmpp 
     821         profdata%nvprot(1)    = ivar1t 
     822         profdata%nvprot(2)    = ivar2t 
     823         profdata%nvprotmpp(1) = ivar1tmpp 
     824         profdata%nvprotmpp(2) = ivar2tmpp 
    800825      ENDIF 
    801826      profdata%nprof        = iprof 
     
    804829      ! Model level search 
    805830      !----------------------------------------------------------------------- 
    806       IF ( ldt3d ) THEN 
     831      IF ( ldvar1 ) THEN 
    807832         CALL obs_level_search( jpk, gdept_1d, & 
    808833            & profdata%nvprot(1), profdata%var(1)%vdep, & 
    809834            & profdata%var(1)%mvk ) 
    810835      ENDIF 
    811       IF ( lds3d ) THEN 
     836      IF ( ldvar2 ) THEN 
    812837         CALL obs_level_search( jpk, gdept_1d, & 
    813838            & profdata%nvprot(2), profdata%var(2)%vdep, & 
    814839            & profdata%var(2)%mvk ) 
    815840      ENDIF 
    816        
     841 
    817842      !----------------------------------------------------------------------- 
    818843      ! Set model equivalent to 99999 
     
    826851      ! Deallocate temporary data 
    827852      !----------------------------------------------------------------------- 
    828       DEALLOCATE( ifileidx, iprofidx, zdat ) 
     853      DEALLOCATE( ifileidx, iprofidx, zdat, clvars ) 
    829854 
    830855      !----------------------------------------------------------------------- 
     
    836861      DEALLOCATE( inpfiles ) 
    837862 
    838    END SUBROUTINE obs_rea_pro_dri 
     863   END SUBROUTINE obs_rea_prof 
    839864 
    840865END MODULE obs_read_prof 
Note: See TracChangeset for help on using the changeset viewer.