- Timestamp:
- 2015-12-04T11:56:46+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_MetOffice_merge_2015/NEMOGCM/NEMO/OPA_SRC/OBS/obs_prep.F90
r5997 r5998 7 7 8 8 !!--------------------------------------------------------------------- 9 !! obs_pre_pro : First level check and screening of T/S profiles 10 !! obs_pre_sla : First level check and screening of SLA observations 11 !! obs_pre_sst : First level check and screening of SLA observations 12 !! obs_pre_seaice : First level check and screening of sea ice observations 13 !! obs_pre_vel : First level check and screening of velocity obs. 14 !! obs_scr : Basic screening of the observations 15 !! obs_coo_tim : Compute number of time steps to the observation time 16 !! obs_sor : Sort the observation arrays 9 !! obs_pre_prof : First level check and screening of profile observations 10 !! obs_pre_surf : First level check and screening of surface observations 11 !! obs_scr : Basic screening of the observations 12 !! obs_coo_tim : Compute number of time steps to the observation time 13 !! obs_sor : Sort the observation arrays 17 14 !!--------------------------------------------------------------------- 18 15 !! * Modules used … … 36 33 37 34 PUBLIC & 38 & obs_pre_pro, & ! First level check and screening of profiles 39 & obs_pre_sla, & ! First level check and screening of SLA data 40 & obs_pre_sst, & ! First level check and screening of SLA data 41 & obs_pre_seaice, & ! First level check and screening of sea ice data 42 & obs_pre_vel, & ! First level check and screening of velocity profiles 43 & calc_month_len ! Calculate the number of days in the months of a year 35 & obs_pre_prof, & ! First level check and screening of profile obs 36 & obs_pre_surf, & ! First level check and screening of surface obs 37 & calc_month_len ! Calculate the number of days in the months of a year 44 38 45 39 !!---------------------------------------------------------------------- … … 49 43 !!---------------------------------------------------------------------- 50 44 51 !! * Substitutions 45 !! * Substitutions 52 46 # include "domzgr_substitute.h90" 53 47 54 48 CONTAINS 55 49 56 SUBROUTINE obs_pre_pro( profdata, prodatqc, ld_t3d, ld_s3d, ld_nea, & 57 & kdailyavtypes ) 58 !!---------------------------------------------------------------------- 59 !! *** ROUTINE obs_pre_pro *** 60 !! 61 !! ** Purpose : First level check and screening of T and S profiles 62 !! 63 !! ** Method : First level check and screening of T and S profiles 64 !! 65 !! ** Action : 66 !! 67 !! References : 68 !! 69 !! History : 70 !! ! 2007-01 (K. Mogensen) Merge of obs_pre_t3d and obs_pre_s3d 71 !! ! 2007-03 (K. Mogensen) General handling of profiles 72 !! ! 2007-06 (K. Mogensen et al) Reject obs. near land. 73 !!---------------------------------------------------------------------- 74 !! * Modules used 75 USE domstp ! Domain: set the time-step 76 USE par_oce ! Ocean parameters 77 USE dom_oce, ONLY : & ! Geographical information 78 & glamt, & 79 & gphit, & 80 & gdept_1d,& 81 & tmask, & 82 & nproc 83 !! * Arguments 84 TYPE(obs_prof), INTENT(INOUT) :: profdata ! Full set of profile data 85 TYPE(obs_prof), INTENT(INOUT) :: prodatqc ! Subset of profile data not failing screening 86 LOGICAL, INTENT(IN) :: ld_t3d ! Switch for temperature 87 LOGICAL, INTENT(IN) :: ld_s3d ! Switch for salinity 88 LOGICAL, INTENT(IN) :: ld_nea ! Switch for rejecting observation near land 89 INTEGER, DIMENSION(imaxavtypes), OPTIONAL :: & 90 & kdailyavtypes! Types for daily averages 91 !! * Local declarations 92 INTEGER :: iyea0 ! Initial date 93 INTEGER :: imon0 ! - (year, month, day, hour, minute) 94 INTEGER :: iday0 95 INTEGER :: ihou0 96 INTEGER :: imin0 97 INTEGER :: icycle ! Current assimilation cycle 98 ! Counters for observations that 99 INTEGER :: iotdobs ! - outside time domain 100 INTEGER :: iosdtobs ! - outside space domain (temperature) 101 INTEGER :: iosdsobs ! - outside space domain (salinity) 102 INTEGER :: ilantobs ! - within a model land cell (temperature) 103 INTEGER :: ilansobs ! - within a model land cell (salinity) 104 INTEGER :: inlatobs ! - close to land (temperature) 105 INTEGER :: inlasobs ! - close to land (salinity) 106 INTEGER :: igrdobs ! - fail the grid search 107 ! Global counters for observations that 108 INTEGER :: iotdobsmpp ! - outside time domain 109 INTEGER :: iosdtobsmpp ! - outside space domain (temperature) 110 INTEGER :: iosdsobsmpp ! - outside space domain (salinity) 111 INTEGER :: ilantobsmpp ! - within a model land cell (temperature) 112 INTEGER :: ilansobsmpp ! - within a model land cell (salinity) 113 INTEGER :: inlatobsmpp ! - close to land (temperature) 114 INTEGER :: inlasobsmpp ! - close to land (salinity) 115 INTEGER :: igrdobsmpp ! - fail the grid search 116 TYPE(obs_prof_valid) :: llvalid ! Profile selection 117 TYPE(obs_prof_valid), DIMENSION(profdata%nvar) :: & 118 & llvvalid ! T,S selection 119 INTEGER :: jvar ! Variable loop variable 120 INTEGER :: jobs ! Obs. loop variable 121 INTEGER :: jstp ! Time loop variable 122 INTEGER :: inrc ! Time index variable 123 124 IF(lwp) WRITE(numout,*)'obs_pre_pro : Preparing the profile observations...' 125 126 ! Initial date initialization (year, month, day, hour, minute) 127 iyea0 = ndate0 / 10000 128 imon0 = ( ndate0 - iyea0 * 10000 ) / 100 129 iday0 = ndate0 - iyea0 * 10000 - imon0 * 100 130 ihou0 = nn_time0 / 100 131 imin0 = ( nn_time0 - ihou0 * 100 ) 132 133 icycle = no ! Assimilation cycle 134 135 ! Diagnotics counters for various failures. 136 137 iotdobs = 0 138 igrdobs = 0 139 iosdtobs = 0 140 iosdsobs = 0 141 ilantobs = 0 142 ilansobs = 0 143 inlatobs = 0 144 inlasobs = 0 145 146 ! ----------------------------------------------------------------------- 147 ! Find time coordinate for profiles 148 ! ----------------------------------------------------------------------- 149 150 IF ( PRESENT(kdailyavtypes) ) THEN 151 CALL obs_coo_tim_prof( icycle, & 152 & iyea0, imon0, iday0, ihou0, imin0, & 153 & profdata%nprof, profdata%nyea, profdata%nmon, & 154 & profdata%nday, profdata%nhou, profdata%nmin, & 155 & profdata%ntyp, profdata%nqc, profdata%mstp, & 156 & iotdobs, kdailyavtypes = kdailyavtypes ) 157 ELSE 158 CALL obs_coo_tim_prof( icycle, & 159 & iyea0, imon0, iday0, ihou0, imin0, & 160 & profdata%nprof, profdata%nyea, profdata%nmon, & 161 & profdata%nday, profdata%nhou, profdata%nmin, & 162 & profdata%ntyp, profdata%nqc, profdata%mstp, & 163 & iotdobs ) 164 ENDIF 165 CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp ) 166 167 ! ----------------------------------------------------------------------- 168 ! Check for profiles failing the grid search 169 ! ----------------------------------------------------------------------- 170 171 CALL obs_coo_grd( profdata%nprof, profdata%mi, profdata%mj, & 172 & profdata%nqc, igrdobs ) 173 174 CALL obs_mpp_sum_integer( igrdobs, igrdobsmpp ) 175 176 ! ----------------------------------------------------------------------- 177 ! Reject all observations for profiles with nqc > 10 178 ! ----------------------------------------------------------------------- 179 180 CALL obs_pro_rej( profdata ) 181 182 ! ----------------------------------------------------------------------- 183 ! Check for land points. This includes points below the model 184 ! bathymetry so this is done for every point in the profile 185 ! ----------------------------------------------------------------------- 186 187 ! Temperature 188 189 CALL obs_coo_spc_3d( profdata%nprof, profdata%nvprot(1), & 190 & profdata%npvsta(:,1), profdata%npvend(:,1), & 191 & jpi, jpj, & 192 & jpk, & 193 & profdata%mi, profdata%mj, & 194 & profdata%var(1)%mvk, & 195 & profdata%rlam, profdata%rphi, & 196 & profdata%var(1)%vdep, & 197 & glamt, gphit, & 198 & gdept_1d, tmask, & 199 & profdata%nqc, profdata%var(1)%nvqc, & 200 & iosdtobs, ilantobs, & 201 & inlatobs, ld_nea ) 202 203 CALL obs_mpp_sum_integer( iosdtobs, iosdtobsmpp ) 204 CALL obs_mpp_sum_integer( ilantobs, ilantobsmpp ) 205 CALL obs_mpp_sum_integer( inlatobs, inlatobsmpp ) 206 207 ! Salinity 208 209 CALL obs_coo_spc_3d( profdata%nprof, profdata%nvprot(2), & 210 & profdata%npvsta(:,2), profdata%npvend(:,2), & 211 & jpi, jpj, & 212 & jpk, & 213 & profdata%mi, profdata%mj, & 214 & profdata%var(2)%mvk, & 215 & profdata%rlam, profdata%rphi, & 216 & profdata%var(2)%vdep, & 217 & glamt, gphit, & 218 & gdept_1d, tmask, & 219 & profdata%nqc, profdata%var(2)%nvqc, & 220 & iosdsobs, ilansobs, & 221 & inlasobs, ld_nea ) 222 223 CALL obs_mpp_sum_integer( iosdsobs, iosdsobsmpp ) 224 CALL obs_mpp_sum_integer( ilansobs, ilansobsmpp ) 225 CALL obs_mpp_sum_integer( inlasobs, inlasobsmpp ) 226 227 ! ----------------------------------------------------------------------- 228 ! Copy useful data from the profdata data structure to 229 ! the prodatqc data structure 230 ! ----------------------------------------------------------------------- 231 232 ! Allocate the selection arrays 233 234 ALLOCATE( llvalid%luse(profdata%nprof) ) 235 DO jvar = 1,profdata%nvar 236 ALLOCATE( llvvalid(jvar)%luse(profdata%nvprot(jvar)) ) 237 END DO 238 239 ! We want all data which has qc flags <= 10 240 241 llvalid%luse(:) = ( profdata%nqc(:) <= 10 ) 242 DO jvar = 1,profdata%nvar 243 llvvalid(jvar)%luse(:) = ( profdata%var(jvar)%nvqc(:) <= 10 ) 244 END DO 245 246 ! The actual copying 247 248 CALL obs_prof_compress( profdata, prodatqc, .TRUE., numout, & 249 & lvalid=llvalid, lvvalid=llvvalid ) 250 251 ! Dellocate the selection arrays 252 DEALLOCATE( llvalid%luse ) 253 DO jvar = 1,profdata%nvar 254 DEALLOCATE( llvvalid(jvar)%luse ) 255 END DO 256 257 ! ----------------------------------------------------------------------- 258 ! Print information about what observations are left after qc 259 ! ----------------------------------------------------------------------- 260 261 ! Update the total observation counter array 262 263 IF(lwp) THEN 264 WRITE(numout,*) 265 WRITE(numout,*) 'obs_pre_pro :' 266 WRITE(numout,*) '~~~~~~~~~~~' 267 WRITE(numout,*) 268 WRITE(numout,*) ' Profiles outside time domain = ', & 269 & iotdobsmpp 270 WRITE(numout,*) ' Remaining profiles that failed grid search = ', & 271 & igrdobsmpp 272 WRITE(numout,*) ' Remaining T data outside space domain = ', & 273 & iosdtobsmpp 274 WRITE(numout,*) ' Remaining T data at land points = ', & 275 & ilantobsmpp 276 IF (ld_nea) THEN 277 WRITE(numout,*) ' Remaining T data near land points (removed) = ',& 278 & inlatobsmpp 279 ELSE 280 WRITE(numout,*) ' Remaining T data near land points (kept) = ',& 281 & inlatobsmpp 282 ENDIF 283 WRITE(numout,*) ' T data accepted = ', & 284 & prodatqc%nvprotmpp(1) 285 WRITE(numout,*) ' Remaining S data outside space domain = ', & 286 & iosdsobsmpp 287 WRITE(numout,*) ' Remaining S data at land points = ', & 288 & ilansobsmpp 289 IF (ld_nea) THEN 290 WRITE(numout,*) ' Remaining S data near land points (removed) = ',& 291 & inlasobsmpp 292 ELSE 293 WRITE(numout,*) ' Remaining S data near land points (kept) = ',& 294 & inlasobsmpp 295 ENDIF 296 WRITE(numout,*) ' S data accepted = ', & 297 & prodatqc%nvprotmpp(2) 298 299 WRITE(numout,*) 300 WRITE(numout,*) ' Number of observations per time step :' 301 WRITE(numout,*) 302 WRITE(numout,997) 303 WRITE(numout,998) 304 ENDIF 305 306 DO jobs = 1, prodatqc%nprof 307 inrc = prodatqc%mstp(jobs) + 2 - nit000 308 prodatqc%npstp(inrc) = prodatqc%npstp(inrc) + 1 309 DO jvar = 1, prodatqc%nvar 310 IF ( prodatqc%npvend(jobs,jvar) > 0 ) THEN 311 prodatqc%nvstp(inrc,jvar) = prodatqc%nvstp(inrc,jvar) + & 312 & ( prodatqc%npvend(jobs,jvar) - & 313 & prodatqc%npvsta(jobs,jvar) + 1 ) 314 ENDIF 315 END DO 316 END DO 317 318 319 CALL obs_mpp_sum_integers( prodatqc%npstp, prodatqc%npstpmpp, & 320 & nitend - nit000 + 2 ) 321 DO jvar = 1, prodatqc%nvar 322 CALL obs_mpp_sum_integers( prodatqc%nvstp(:,jvar), & 323 & prodatqc%nvstpmpp(:,jvar), & 324 & nitend - nit000 + 2 ) 325 END DO 326 327 IF ( lwp ) THEN 328 DO jstp = nit000 - 1, nitend 329 inrc = jstp - nit000 + 2 330 WRITE(numout,999) jstp, prodatqc%npstpmpp(inrc), & 331 & prodatqc%nvstpmpp(inrc,1), & 332 & prodatqc%nvstpmpp(inrc,2) 333 END DO 334 ENDIF 335 336 997 FORMAT(10X,'Time step',5X,'Profiles',5X,'Temperature',5X,'Salinity') 337 998 FORMAT(10X,'---------',5X,'--------',5X,'-----------',5X,'--------') 338 999 FORMAT(10X,I9,5X,I8,5X,I11,5X,I8) 339 340 END SUBROUTINE obs_pre_pro 341 342 SUBROUTINE obs_pre_sla( sladata, sladatqc, ld_sla, ld_nea ) 50 SUBROUTINE obs_pre_surf( surfdata, surfdataqc, ld_nea ) 343 51 !!---------------------------------------------------------------------- 344 52 !! *** ROUTINE obs_pre_sla *** 345 53 !! 346 !! ** Purpose : First level check and screening of SLAobservations347 !! 348 !! ** Method : First level check and screening of SLAobservations54 !! ** Purpose : First level check and screening of surface observations 55 !! 56 !! ** Method : First level check and screening of surface observations 349 57 !! 350 58 !! ** Action : … … 355 63 !! ! 2007-03 (A. Weaver, K. Mogensen) Original 356 64 !! ! 2007-06 (K. Mogensen et al) Reject obs. near land. 65 !! ! 2015-02 (M. Martin) Combined routine for surface types. 357 66 !!---------------------------------------------------------------------- 358 67 !! * Modules used … … 365 74 & nproc 366 75 !! * Arguments 367 TYPE(obs_surf), INTENT(INOUT) :: sladata ! Full set of SLA data 368 TYPE(obs_surf), INTENT(INOUT) :: sladatqc ! Subset of SLA data not failing screening 369 LOGICAL, INTENT(IN) :: ld_sla ! Switch for SLA data 76 TYPE(obs_surf), INTENT(INOUT) :: surfdata ! Full set of surface data 77 TYPE(obs_surf), INTENT(INOUT) :: surfdataqc ! Subset of surface data not failing screening 370 78 LOGICAL, INTENT(IN) :: ld_nea ! Switch for rejecting observation near land 371 79 !! * Local declarations … … 394 102 INTEGER :: inrc ! Time index variable 395 103 396 IF(lwp) WRITE(numout,*)'obs_pre_sla : Preparing the SLA observations...' 397 104 IF(lwp) WRITE(numout,*)'obs_pre_surf : Preparing the surface observations...' 105 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 106 398 107 ! Initial date initialization (year, month, day, hour, minute) 399 108 iyea0 = ndate0 / 10000 400 109 imon0 = ( ndate0 - iyea0 * 10000 ) / 100 401 110 iday0 = ndate0 - iyea0 * 10000 - imon0 * 100 402 ihou0 = 111 ihou0 = nn_time0 / 100 403 112 imin0 = ( nn_time0 - ihou0 * 100 ) 404 113 … … 414 123 415 124 ! ----------------------------------------------------------------------- 416 ! Find time coordinate for SLAdata125 ! Find time coordinate for surface data 417 126 ! ----------------------------------------------------------------------- 418 127 419 128 CALL obs_coo_tim( icycle, & 420 129 & iyea0, imon0, iday0, ihou0, imin0, & 421 & s ladata%nsurf, sladata%nyea, sladata%nmon, &422 & s ladata%nday, sladata%nhou, sladata%nmin, &423 & s ladata%nqc, sladata%mstp, iotdobs )130 & surfdata%nsurf, surfdata%nyea, surfdata%nmon, & 131 & surfdata%nday, surfdata%nhou, surfdata%nmin, & 132 & surfdata%nqc, surfdata%mstp, iotdobs ) 424 133 425 134 CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp ) 426 135 427 136 ! ----------------------------------------------------------------------- 428 ! Check for SLAdata failing the grid search429 ! ----------------------------------------------------------------------- 430 431 CALL obs_coo_grd( s ladata%nsurf, sladata%mi, sladata%mj, &432 & s ladata%nqc, igrdobs )137 ! Check for surface data failing the grid search 138 ! ----------------------------------------------------------------------- 139 140 CALL obs_coo_grd( surfdata%nsurf, surfdata%mi, surfdata%mj, & 141 & surfdata%nqc, igrdobs ) 433 142 434 143 CALL obs_mpp_sum_integer( igrdobs, igrdobsmpp ) … … 438 147 ! ----------------------------------------------------------------------- 439 148 440 CALL obs_coo_spc_2d( s ladata%nsurf, &149 CALL obs_coo_spc_2d( surfdata%nsurf, & 441 150 & jpi, jpj, & 442 & s ladata%mi, sladata%mj, &443 & s ladata%rlam, sladata%rphi, &151 & surfdata%mi, surfdata%mj, & 152 & surfdata%rlam, surfdata%rphi, & 444 153 & glamt, gphit, & 445 & tmask(:,:,1), s ladata%nqc, &154 & tmask(:,:,1), surfdata%nqc, & 446 155 & iosdsobs, ilansobs, & 447 156 & inlasobs, ld_nea ) … … 452 161 453 162 ! ----------------------------------------------------------------------- 454 ! Copy useful data from the s ladata data structure to455 ! the s ladatqc data structure163 ! Copy useful data from the surfdata data structure to 164 ! the surfdataqc data structure 456 165 ! ----------------------------------------------------------------------- 457 166 458 167 ! Allocate the selection arrays 459 168 460 ALLOCATE( llvalid(s ladata%nsurf) )169 ALLOCATE( llvalid(surfdata%nsurf) ) 461 170 462 171 ! We want all data which has qc flags <= 10 463 172 464 llvalid(:) = ( s ladata%nqc(:) <= 10 )173 llvalid(:) = ( surfdata%nqc(:) <= 10 ) 465 174 466 175 ! The actual copying 467 176 468 CALL obs_surf_compress( s ladata, sladatqc, .TRUE., numout, &177 CALL obs_surf_compress( surfdata, surfdataqc, .TRUE., numout, & 469 178 & lvalid=llvalid ) 470 179 … … 480 189 IF(lwp) THEN 481 190 WRITE(numout,*) 482 WRITE(numout,*) 'obs_pre_sla :' 483 WRITE(numout,*) '~~~~~~~~~~~' 484 WRITE(numout,*) 485 WRITE(numout,*) ' SLA data outside time domain = ', & 191 WRITE(numout,*) ' '//surfdataqc%cvars(1)//' data outside time domain = ', & 486 192 & iotdobsmpp 487 WRITE(numout,*) ' Remaining SLAdata that failed grid search = ', &193 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data that failed grid search = ', & 488 194 & igrdobsmpp 489 WRITE(numout,*) ' Remaining SLAdata outside space domain = ', &195 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data outside space domain = ', & 490 196 & iosdsobsmpp 491 WRITE(numout,*) ' Remaining SLAdata at land points = ', &197 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data at land points = ', & 492 198 & ilansobsmpp 493 199 IF (ld_nea) THEN 494 WRITE(numout,*) ' Remaining SLAdata near land points (removed) = ', &200 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data near land points (removed) = ', & 495 201 & inlasobsmpp 496 202 ELSE 497 WRITE(numout,*) ' Remaining SLAdata near land points (kept) = ', &203 WRITE(numout,*) ' Remaining '//surfdataqc%cvars(1)//' data near land points (kept) = ', & 498 204 & inlasobsmpp 499 205 ENDIF 500 WRITE(numout,*) ' SLAdata accepted = ', &501 & s ladatqc%nsurfmpp206 WRITE(numout,*) ' '//surfdataqc%cvars(1)//' data accepted = ', & 207 & surfdataqc%nsurfmpp 502 208 503 209 WRITE(numout,*) 504 210 WRITE(numout,*) ' Number of observations per time step :' 505 211 WRITE(numout,*) 506 WRITE(numout,1997) 507 WRITE(numout,1998) 212 WRITE(numout,'(10X,A,10X,A)')'Time step',surfdataqc%cvars(1) 213 WRITE(numout,'(10X,A,5X,A)')'---------','-----------------' 214 CALL FLUSH(numout) 508 215 ENDIF 509 216 510 DO jobs = 1, s ladatqc%nsurf511 inrc = s ladatqc%mstp(jobs) + 2 - nit000512 s ladatqc%nsstp(inrc) = sladatqc%nsstp(inrc) + 1217 DO jobs = 1, surfdataqc%nsurf 218 inrc = surfdataqc%mstp(jobs) + 2 - nit000 219 surfdataqc%nsstp(inrc) = surfdataqc%nsstp(inrc) + 1 513 220 END DO 514 221 515 CALL obs_mpp_sum_integers( s ladatqc%nsstp, sladatqc%nsstpmpp, &222 CALL obs_mpp_sum_integers( surfdataqc%nsstp, surfdataqc%nsstpmpp, & 516 223 & nitend - nit000 + 2 ) 517 224 … … 519 226 DO jstp = nit000 - 1, nitend 520 227 inrc = jstp - nit000 + 2 521 WRITE(numout,1999) jstp, sladatqc%nsstpmpp(inrc) 228 WRITE(numout,1999) jstp, surfdataqc%nsstpmpp(inrc) 229 CALL FLUSH(numout) 522 230 END DO 523 231 ENDIF 524 232 525 1997 FORMAT(10X,'Time step',5X,'Sea level anomaly')526 1998 FORMAT(10X,'---------',5X,'-----------------')527 233 1999 FORMAT(10X,I9,5X,I17) 528 234 529 END SUBROUTINE obs_pre_sla 530 531 SUBROUTINE obs_pre_sst( sstdata, sstdatqc, ld_sst, ld_nea ) 532 !!---------------------------------------------------------------------- 533 !! *** ROUTINE obs_pre_sst *** 534 !! 535 !! ** Purpose : First level check and screening of SST observations 536 !! 537 !! ** Method : First level check and screening of SST observations 538 !! 539 !! ** Action : 540 !! 541 !! References : 542 !! 235 END SUBROUTINE obs_pre_surf 236 237 238 SUBROUTINE obs_pre_prof( profdata, prodatqc, ld_var1, ld_var2, & 239 & kpi, kpj, kpk, & 240 & zmask1, pglam1, pgphi1, zmask2, pglam2, pgphi2, & 241 & ld_nea, kdailyavtypes ) 242 243 !!---------------------------------------------------------------------- 244 !! *** ROUTINE obs_pre_prof *** 245 !! 246 !! ** Purpose : First level check and screening of profiles 247 !! 248 !! ** Method : First level check and screening of profiles 249 !! 543 250 !! History : 544 !! ! 2007-03 (S. Ricci) SST data preparation 251 !! ! 2007-06 (K. Mogensen) original : T and S profile data 252 !! ! 2008-09 (M. Valdivieso) : TAO velocity data 253 !! ! 2009-01 (K. Mogensen) : New feedback stricture 254 !! ! 2015-02 (M. Martin) : Combined profile routine. 255 !! 545 256 !!---------------------------------------------------------------------- 546 257 !! * Modules used … … 548 259 USE par_oce ! Ocean parameters 549 260 USE dom_oce, ONLY : & ! Geographical information 550 & glamt, & 551 & gphit, & 552 & tmask, & 261 & gdept_1d, & 553 262 & nproc 554 !! * Arguments 555 TYPE(obs_surf), INTENT(INOUT) :: sstdata ! Full set of SST data 556 TYPE(obs_surf), INTENT(INOUT) :: sstdatqc ! Subset of SST data not failing screening 557 LOGICAL :: ld_sst ! Switch for SST data 558 LOGICAL :: ld_nea ! Switch for rejecting observation near land 559 !! * Local declarations 560 INTEGER :: iyea0 ! Initial date 561 INTEGER :: imon0 ! - (year, month, day, hour, minute) 562 INTEGER :: iday0 563 INTEGER :: ihou0 564 INTEGER :: imin0 565 INTEGER :: icycle ! Current assimilation cycle 566 ! Counters for observations that 567 INTEGER :: iotdobs ! - outside time domain 568 INTEGER :: iosdsobs ! - outside space domain 569 INTEGER :: ilansobs ! - within a model land cell 570 INTEGER :: inlasobs ! - close to land 571 INTEGER :: igrdobs ! - fail the grid search 572 ! Global counters for observations that 573 INTEGER :: iotdobsmpp ! - outside time domain 574 INTEGER :: iosdsobsmpp ! - outside space domain 575 INTEGER :: ilansobsmpp ! - within a model land cell 576 INTEGER :: inlasobsmpp ! - close to land 577 INTEGER :: igrdobsmpp ! - fail the grid search 578 LOGICAL, DIMENSION(:), ALLOCATABLE :: & 579 & llvalid ! SST data selection 580 INTEGER :: jobs ! Obs. loop variable 581 INTEGER :: jstp ! Time loop variable 582 INTEGER :: inrc ! Time index variable 583 584 IF(lwp) WRITE(numout,*)'obs_pre_sst : Preparing the SST observations...' 585 586 ! Initial date initialization (year, month, day, hour, minute) 587 iyea0 = ndate0 / 10000 588 imon0 = ( ndate0 - iyea0 * 10000 ) / 100 589 iday0 = ndate0 - iyea0 * 10000 - imon0 * 100 590 ihou0 = nn_time0 / 100 591 imin0 = ( nn_time0 - ihou0 * 100 ) 592 593 icycle = no ! Assimilation cycle 594 595 ! Diagnotics counters for various failures. 596 597 iotdobs = 0 598 igrdobs = 0 599 iosdsobs = 0 600 ilansobs = 0 601 inlasobs = 0 602 603 ! ----------------------------------------------------------------------- 604 ! Find time coordinate for SST data 605 ! ----------------------------------------------------------------------- 606 607 CALL obs_coo_tim( icycle, & 608 & iyea0, imon0, iday0, ihou0, imin0, & 609 & sstdata%nsurf, sstdata%nyea, sstdata%nmon, & 610 & sstdata%nday, sstdata%nhou, sstdata%nmin, & 611 & sstdata%nqc, sstdata%mstp, iotdobs ) 612 CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp ) 613 ! ----------------------------------------------------------------------- 614 ! Check for SST data failing the grid search 615 ! ----------------------------------------------------------------------- 616 617 CALL obs_coo_grd( sstdata%nsurf, sstdata%mi, sstdata%mj, & 618 & sstdata%nqc, igrdobs ) 619 CALL obs_mpp_sum_integer( igrdobs, igrdobsmpp ) 620 621 ! ----------------------------------------------------------------------- 622 ! Check for land points. 623 ! ----------------------------------------------------------------------- 624 625 CALL obs_coo_spc_2d( sstdata%nsurf, & 626 & jpi, jpj, & 627 & sstdata%mi, sstdata%mj, & 628 & sstdata%rlam, sstdata%rphi, & 629 & glamt, gphit, & 630 & tmask(:,:,1), sstdata%nqc, & 631 & iosdsobs, ilansobs, & 632 & inlasobs, ld_nea ) 633 634 CALL obs_mpp_sum_integer( iosdsobs, iosdsobsmpp ) 635 CALL obs_mpp_sum_integer( ilansobs, ilansobsmpp ) 636 CALL obs_mpp_sum_integer( inlasobs, inlasobsmpp ) 637 638 ! ----------------------------------------------------------------------- 639 ! Copy useful data from the sstdata data structure to 640 ! the sstdatqc data structure 641 ! ----------------------------------------------------------------------- 642 643 ! Allocate the selection arrays 644 645 ALLOCATE( llvalid(sstdata%nsurf) ) 646 647 ! We want all data which has qc flags <= 0 648 649 llvalid(:) = ( sstdata%nqc(:) <= 10 ) 650 651 ! The actual copying 652 653 CALL obs_surf_compress( sstdata, sstdatqc, .TRUE., numout, & 654 & lvalid=llvalid ) 655 656 ! Dellocate the selection arrays 657 DEALLOCATE( llvalid ) 658 659 ! ----------------------------------------------------------------------- 660 ! Print information about what observations are left after qc 661 ! ----------------------------------------------------------------------- 662 663 ! Update the total observation counter array 664 665 IF(lwp) THEN 666 WRITE(numout,*) 667 WRITE(numout,*) 'obs_pre_sst :' 668 WRITE(numout,*) '~~~~~~~~~~~' 669 WRITE(numout,*) 670 WRITE(numout,*) ' SST data outside time domain = ', & 671 & iotdobsmpp 672 WRITE(numout,*) ' Remaining SST data that failed grid search = ', & 673 & igrdobsmpp 674 WRITE(numout,*) ' Remaining SST data outside space domain = ', & 675 & iosdsobsmpp 676 WRITE(numout,*) ' Remaining SST data at land points = ', & 677 & ilansobsmpp 678 IF (ld_nea) THEN 679 WRITE(numout,*) ' Remaining SST data near land points (removed) = ', & 680 & inlasobsmpp 681 ELSE 682 WRITE(numout,*) ' Remaining SST data near land points (kept) = ', & 683 & inlasobsmpp 684 ENDIF 685 WRITE(numout,*) ' SST data accepted = ', & 686 & sstdatqc%nsurfmpp 687 688 WRITE(numout,*) 689 WRITE(numout,*) ' Number of observations per time step :' 690 WRITE(numout,*) 691 WRITE(numout,1997) 692 WRITE(numout,1998) 693 ENDIF 694 695 DO jobs = 1, sstdatqc%nsurf 696 inrc = sstdatqc%mstp(jobs) + 2 - nit000 697 sstdatqc%nsstp(inrc) = sstdatqc%nsstp(inrc) + 1 698 END DO 699 700 CALL obs_mpp_sum_integers( sstdatqc%nsstp, sstdatqc%nsstpmpp, & 701 & nitend - nit000 + 2 ) 702 703 IF ( lwp ) THEN 704 DO jstp = nit000 - 1, nitend 705 inrc = jstp - nit000 + 2 706 WRITE(numout,1999) jstp, sstdatqc%nsstpmpp(inrc) 707 END DO 708 ENDIF 709 710 1997 FORMAT(10X,'Time step',5X,'Sea surface temperature') 711 1998 FORMAT(10X,'---------',5X,'-----------------') 712 1999 FORMAT(10X,I9,5X,I17) 713 714 END SUBROUTINE obs_pre_sst 715 716 SUBROUTINE obs_pre_seaice( seaicedata, seaicedatqc, ld_seaice, ld_nea ) 717 !!---------------------------------------------------------------------- 718 !! *** ROUTINE obs_pre_seaice *** 719 !! 720 !! ** Purpose : First level check and screening of Sea Ice observations 721 !! 722 !! ** Method : First level check and screening of Sea Ice observations 723 !! 724 !! ** Action : 725 !! 726 !! References : 727 !! 728 !! History : 729 !! ! 2007-11 (D. Lea) based on obs_pre_sst 730 !!---------------------------------------------------------------------- 731 !! * Modules used 732 USE domstp ! Domain: set the time-step 733 USE par_oce ! Ocean parameters 734 USE dom_oce, ONLY : & ! Geographical information 735 & glamt, & 736 & gphit, & 737 & tmask, & 738 & nproc 739 !! * Arguments 740 TYPE(obs_surf), INTENT(INOUT) :: seaicedata ! Full set of Sea Ice data 741 TYPE(obs_surf), INTENT(INOUT) :: seaicedatqc ! Subset of sea ice data not failing screening 742 LOGICAL :: ld_seaice ! Switch for sea ice data 743 LOGICAL :: ld_nea ! Switch for rejecting observation near land 744 !! * Local declarations 745 INTEGER :: iyea0 ! Initial date 746 INTEGER :: imon0 ! - (year, month, day, hour, minute) 747 INTEGER :: iday0 748 INTEGER :: ihou0 749 INTEGER :: imin0 750 INTEGER :: icycle ! Current assimilation cycle 751 ! Counters for observations that 752 INTEGER :: iotdobs ! - outside time domain 753 INTEGER :: iosdsobs ! - outside space domain 754 INTEGER :: ilansobs ! - within a model land cell 755 INTEGER :: inlasobs ! - close to land 756 INTEGER :: igrdobs ! - fail the grid search 757 ! Global counters for observations that 758 INTEGER :: iotdobsmpp ! - outside time domain 759 INTEGER :: iosdsobsmpp ! - outside space domain 760 INTEGER :: ilansobsmpp ! - within a model land cell 761 INTEGER :: inlasobsmpp ! - close to land 762 INTEGER :: igrdobsmpp ! - fail the grid search 763 LOGICAL, DIMENSION(:), ALLOCATABLE :: & 764 & llvalid ! data selection 765 INTEGER :: jobs ! Obs. loop variable 766 INTEGER :: jstp ! Time loop variable 767 INTEGER :: inrc ! Time index variable 768 769 IF (lwp) WRITE(numout,*)'obs_pre_seaice : Preparing the sea ice observations...' 770 771 ! Initial date initialization (year, month, day, hour, minute) 772 iyea0 = ndate0 / 10000 773 imon0 = ( ndate0 - iyea0 * 10000 ) / 100 774 iday0 = ndate0 - iyea0 * 10000 - imon0 * 100 775 ihou0 = nn_time0 / 100 776 imin0 = ( nn_time0 - ihou0 * 100 ) 777 778 icycle = no ! Assimilation cycle 779 780 ! Diagnotics counters for various failures. 781 782 iotdobs = 0 783 igrdobs = 0 784 iosdsobs = 0 785 ilansobs = 0 786 inlasobs = 0 787 788 ! ----------------------------------------------------------------------- 789 ! Find time coordinate for sea ice data 790 ! ----------------------------------------------------------------------- 791 792 CALL obs_coo_tim( icycle, & 793 & iyea0, imon0, iday0, ihou0, imin0, & 794 & seaicedata%nsurf, seaicedata%nyea, seaicedata%nmon, & 795 & seaicedata%nday, seaicedata%nhou, seaicedata%nmin, & 796 & seaicedata%nqc, seaicedata%mstp, iotdobs ) 797 CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp ) 798 ! ----------------------------------------------------------------------- 799 ! Check for sea ice data failing the grid search 800 ! ----------------------------------------------------------------------- 801 802 CALL obs_coo_grd( seaicedata%nsurf, seaicedata%mi, seaicedata%mj, & 803 & seaicedata%nqc, igrdobs ) 804 CALL obs_mpp_sum_integer( igrdobs, igrdobsmpp ) 805 806 ! ----------------------------------------------------------------------- 807 ! Check for land points. 808 ! ----------------------------------------------------------------------- 809 810 CALL obs_coo_spc_2d( seaicedata%nsurf, & 811 & jpi, jpj, & 812 & seaicedata%mi, seaicedata%mj, & 813 & seaicedata%rlam, seaicedata%rphi, & 814 & glamt, gphit, & 815 & tmask(:,:,1), seaicedata%nqc, & 816 & iosdsobs, ilansobs, & 817 & inlasobs, ld_nea ) 818 819 CALL obs_mpp_sum_integer( iosdsobs, iosdsobsmpp ) 820 CALL obs_mpp_sum_integer( ilansobs, ilansobsmpp ) 821 CALL obs_mpp_sum_integer( inlasobs, inlasobsmpp ) 822 823 ! ----------------------------------------------------------------------- 824 ! Copy useful data from the seaicedata data structure to 825 ! the seaicedatqc data structure 826 ! ----------------------------------------------------------------------- 827 828 ! Allocate the selection arrays 829 830 ALLOCATE( llvalid(seaicedata%nsurf) ) 831 832 ! We want all data which has qc flags <= 0 833 834 llvalid(:) = ( seaicedata%nqc(:) <= 10 ) 835 836 ! The actual copying 837 838 CALL obs_surf_compress( seaicedata, seaicedatqc, .TRUE., numout, & 839 & lvalid=llvalid ) 840 841 ! Dellocate the selection arrays 842 DEALLOCATE( llvalid ) 843 844 ! ----------------------------------------------------------------------- 845 ! Print information about what observations are left after qc 846 ! ----------------------------------------------------------------------- 847 848 ! Update the total observation counter array 849 850 IF(lwp) THEN 851 WRITE(numout,*) 852 WRITE(numout,*) 'obs_pre_seaice :' 853 WRITE(numout,*) '~~~~~~~~~~~' 854 WRITE(numout,*) 855 WRITE(numout,*) ' Sea ice data outside time domain = ', & 856 & iotdobsmpp 857 WRITE(numout,*) ' Remaining sea ice data that failed grid search = ', & 858 & igrdobsmpp 859 WRITE(numout,*) ' Remaining sea ice data outside space domain = ', & 860 & iosdsobsmpp 861 WRITE(numout,*) ' Remaining sea ice data at land points = ', & 862 & ilansobsmpp 863 IF (ld_nea) THEN 864 WRITE(numout,*) ' Remaining sea ice data near land points (removed) = ', & 865 & inlasobsmpp 866 ELSE 867 WRITE(numout,*) ' Remaining sea ice data near land points (kept) = ', & 868 & inlasobsmpp 869 ENDIF 870 WRITE(numout,*) ' Sea ice data accepted = ', & 871 & seaicedatqc%nsurfmpp 872 873 WRITE(numout,*) 874 WRITE(numout,*) ' Number of observations per time step :' 875 WRITE(numout,*) 876 WRITE(numout,1997) 877 WRITE(numout,1998) 878 ENDIF 879 880 DO jobs = 1, seaicedatqc%nsurf 881 inrc = seaicedatqc%mstp(jobs) + 2 - nit000 882 seaicedatqc%nsstp(inrc) = seaicedatqc%nsstp(inrc) + 1 883 END DO 884 885 CALL obs_mpp_sum_integers( seaicedatqc%nsstp, seaicedatqc%nsstpmpp, & 886 & nitend - nit000 + 2 ) 887 888 IF ( lwp ) THEN 889 DO jstp = nit000 - 1, nitend 890 inrc = jstp - nit000 + 2 891 WRITE(numout,1999) jstp, seaicedatqc%nsstpmpp(inrc) 892 END DO 893 ENDIF 894 895 1997 FORMAT(10X,'Time step',5X,'Sea ice data ') 896 1998 FORMAT(10X,'---------',5X,'-----------------') 897 1999 FORMAT(10X,I9,5X,I17) 898 899 END SUBROUTINE obs_pre_seaice 900 901 SUBROUTINE obs_pre_vel( profdata, prodatqc, ld_vel3d, ld_nea, ld_dailyav ) 902 !!---------------------------------------------------------------------- 903 !! *** ROUTINE obs_pre_taovel *** 904 !! 905 !! ** Purpose : First level check and screening of U and V profiles 906 !! 907 !! ** Method : First level check and screening of U and V profiles 908 !! 909 !! History : 910 !! ! 2007-06 (K. Mogensen) original : T and S profile data 911 !! ! 2008-09 (M. Valdivieso) : TAO velocity data 912 !! ! 2009-01 (K. Mogensen) : New feedback strictuer 913 !! 914 !!---------------------------------------------------------------------- 915 !! * Modules used 916 USE domstp ! Domain: set the time-step 917 USE par_oce ! Ocean parameters 918 USE dom_oce, ONLY : & ! Geographical information 919 & glamt, glamu, glamv, & 920 & gphit, gphiu, gphiv, & 921 & gdept_1d, & 922 & tmask, umask, vmask, & 923 & nproc 263 924 264 !! * Arguments 925 265 TYPE(obs_prof), INTENT(INOUT) :: profdata ! Full set of profile data 926 266 TYPE(obs_prof), INTENT(INOUT) :: prodatqc ! Subset of profile data not failing screening 927 LOGICAL, INTENT(IN) :: ld_vel3d ! Switch for zonal and meridional velocity components 928 LOGICAL, INTENT(IN) :: ld_nea ! Switch for rejecting observation near land 929 LOGICAL, INTENT(IN) :: ld_dailyav ! Switch for daily average data 267 LOGICAL, INTENT(IN) :: ld_var1 ! Observed variables switches 268 LOGICAL, INTENT(IN) :: ld_var2 269 LOGICAL, INTENT(IN) :: ld_nea ! Switch for rejecting observation near land 270 INTEGER, INTENT(IN) :: kpi, kpj, kpk ! Local domain sizes 271 INTEGER, DIMENSION(imaxavtypes), OPTIONAL :: & 272 & kdailyavtypes ! Types for daily averages 273 REAL(wp), INTENT(IN), DIMENSION(kpi,kpj,kpk) :: & 274 & zmask1, & 275 & zmask2 276 REAL(wp), INTENT(IN), DIMENSION(kpi,kpj) :: & 277 & pglam1, & 278 & pglam2, & 279 & pgphi1, & 280 & pgphi2 281 930 282 !! * Local declarations 931 283 INTEGER :: iyea0 ! Initial date … … 935 287 INTEGER :: imin0 936 288 INTEGER :: icycle ! Current assimilation cycle 937 ! Counters for observations that 289 ! Counters for observations that are 938 290 INTEGER :: iotdobs ! - outside time domain 939 INTEGER :: iosd uobs ! - outside space domain (zonal velocity component)940 INTEGER :: iosdv obs ! - outside space domain (meridional velocity component)941 INTEGER :: ilan uobs ! - within a model land cell (zonal velocity component)942 INTEGER :: ilanv obs ! - within a model land cell (meridional velocity component)943 INTEGER :: inla uobs ! - close to land (zonal velocity component)944 INTEGER :: inlav obs ! - close to land (meridional velocity component)291 INTEGER :: iosdv1obs ! - outside space domain (variable 1) 292 INTEGER :: iosdv2obs ! - outside space domain (variable 2) 293 INTEGER :: ilanv1obs ! - within a model land cell (variable 1) 294 INTEGER :: ilanv2obs ! - within a model land cell (variable 2) 295 INTEGER :: inlav1obs ! - close to land (variable 1) 296 INTEGER :: inlav2obs ! - close to land (variable 2) 945 297 INTEGER :: igrdobs ! - fail the grid search 946 298 INTEGER :: iuvchku ! - reject u if v rejected and vice versa 947 299 INTEGER :: iuvchkv ! 948 ! Global counters for observations that 300 ! Global counters for observations that are 949 301 INTEGER :: iotdobsmpp ! - outside time domain 950 INTEGER :: iosd uobsmpp ! - outside space domain (zonal velocity component)951 INTEGER :: iosdv obsmpp ! - outside space domain (meridional velocity component)952 INTEGER :: ilan uobsmpp ! - within a model land cell (zonal velocity component)953 INTEGER :: ilanv obsmpp ! - within a model land cell (meridional velocity component)954 INTEGER :: inla uobsmpp ! - close to land (zonal velocity component)955 INTEGER :: inlav obsmpp ! - close to land (meridional velocity component)302 INTEGER :: iosdv1obsmpp ! - outside space domain (variable 1) 303 INTEGER :: iosdv2obsmpp ! - outside space domain (variable 2) 304 INTEGER :: ilanv1obsmpp ! - within a model land cell (variable 1) 305 INTEGER :: ilanv2obsmpp ! - within a model land cell (variable 2) 306 INTEGER :: inlav1obsmpp ! - close to land (variable 1) 307 INTEGER :: inlav2obsmpp ! - close to land (variable 2) 956 308 INTEGER :: igrdobsmpp ! - fail the grid search 957 INTEGER :: iuvchkumpp ! - reject u if vrejected and vice versa309 INTEGER :: iuvchkumpp ! - reject var1 if var2 rejected and vice versa 958 310 INTEGER :: iuvchkvmpp ! 959 311 TYPE(obs_prof_valid) :: llvalid ! Profile selection 960 312 TYPE(obs_prof_valid), DIMENSION(profdata%nvar) :: & 961 & llvvalid ! U,Vselection313 & llvvalid ! var1,var2 selection 962 314 INTEGER :: jvar ! Variable loop variable 963 315 INTEGER :: jobs ! Obs. loop variable … … 965 317 INTEGER :: inrc ! Time index variable 966 318 967 IF(lwp) WRITE(numout,*)'obs_pre_vel: Preparing the velocity profile data' 319 IF(lwp) WRITE(numout,*)'obs_pre_prof: Preparing the profile data...' 320 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 968 321 969 322 ! Initial date initialization (year, month, day, hour, minute) … … 971 324 imon0 = ( ndate0 - iyea0 * 10000 ) / 100 972 325 iday0 = ndate0 - iyea0 * 10000 - imon0 * 100 973 ihou0 = 326 ihou0 = nn_time0 / 100 974 327 imin0 = ( nn_time0 - ihou0 * 100 ) 975 328 … … 980 333 iotdobs = 0 981 334 igrdobs = 0 982 iosd uobs = 0983 iosdv obs = 0984 ilan uobs = 0985 ilanv obs = 0986 inla uobs = 0987 inlav obs = 0335 iosdv1obs = 0 336 iosdv2obs = 0 337 ilanv1obs = 0 338 ilanv2obs = 0 339 inlav1obs = 0 340 inlav2obs = 0 988 341 iuvchku = 0 989 342 iuvchkv = 0 … … 993 346 ! ----------------------------------------------------------------------- 994 347 995 CALL obs_coo_tim_prof( icycle, & 996 & iyea0, imon0, iday0, ihou0, imin0, & 997 & profdata%nprof, profdata%nyea, profdata%nmon, & 998 & profdata%nday, profdata%nhou, profdata%nmin, & 999 & profdata%ntyp, profdata%nqc, profdata%mstp, & 1000 & iotdobs, ld_dailyav = ld_dailyav ) 1001 348 IF ( PRESENT(kdailyavtypes) ) THEN 349 CALL obs_coo_tim_prof( icycle, & 350 & iyea0, imon0, iday0, ihou0, imin0, & 351 & profdata%nprof, profdata%nyea, profdata%nmon, & 352 & profdata%nday, profdata%nhou, profdata%nmin, & 353 & profdata%ntyp, profdata%nqc, profdata%mstp, & 354 & iotdobs, kdailyavtypes = kdailyavtypes ) 355 ELSE 356 CALL obs_coo_tim_prof( icycle, & 357 & iyea0, imon0, iday0, ihou0, imin0, & 358 & profdata%nprof, profdata%nyea, profdata%nmon, & 359 & profdata%nday, profdata%nhou, profdata%nmin, & 360 & profdata%ntyp, profdata%nqc, profdata%mstp, & 361 & iotdobs ) 362 ENDIF 363 1002 364 CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp ) 1003 365 … … 1024 386 ! ----------------------------------------------------------------------- 1025 387 1026 ! Zonal Velocity Component 1027 388 ! Variable 1 1028 389 CALL obs_coo_spc_3d( profdata%nprof, profdata%nvprot(1), & 1029 390 & profdata%npvsta(:,1), profdata%npvend(:,1), & 1030 391 & jpi, jpj, & 1031 392 & jpk, & 1032 & profdata%mi, profdata%mj, & 393 & profdata%mi, profdata%mj, & 1033 394 & profdata%var(1)%mvk, & 1034 395 & profdata%rlam, profdata%rphi, & 1035 396 & profdata%var(1)%vdep, & 1036 & glamu, gphiu,&1037 & gdept_1d, umask,&397 & pglam1, pgphi1, & 398 & gdept_1d, zmask1, & 1038 399 & profdata%nqc, profdata%var(1)%nvqc, & 1039 & iosduobs, ilanuobs, & 1040 & inlauobs, ld_nea ) 1041 1042 CALL obs_mpp_sum_integer( iosduobs, iosduobsmpp ) 1043 CALL obs_mpp_sum_integer( ilanuobs, ilanuobsmpp ) 1044 CALL obs_mpp_sum_integer( inlauobs, inlauobsmpp ) 1045 1046 ! Meridional Velocity Component 1047 400 & iosdv1obs, ilanv1obs, & 401 & inlav1obs, ld_nea ) 402 403 CALL obs_mpp_sum_integer( iosdv1obs, iosdv1obsmpp ) 404 CALL obs_mpp_sum_integer( ilanv1obs, ilanv1obsmpp ) 405 CALL obs_mpp_sum_integer( inlav1obs, inlav1obsmpp ) 406 407 ! Variable 2 1048 408 CALL obs_coo_spc_3d( profdata%nprof, profdata%nvprot(2), & 1049 409 & profdata%npvsta(:,2), profdata%npvend(:,2), & … … 1054 414 & profdata%rlam, profdata%rphi, & 1055 415 & profdata%var(2)%vdep, & 1056 & glamv, gphiv,&1057 & gdept_1d, vmask,&416 & pglam2, pgphi2, & 417 & gdept_1d, zmask2, & 1058 418 & profdata%nqc, profdata%var(2)%nvqc, & 1059 & iosdv obs, ilanvobs,&1060 & inlav obs, ld_nea )1061 1062 CALL obs_mpp_sum_integer( iosdv obs, iosdvobsmpp )1063 CALL obs_mpp_sum_integer( ilanv obs, ilanvobsmpp )1064 CALL obs_mpp_sum_integer( inlav obs, inlavobsmpp )419 & iosdv2obs, ilanv2obs, & 420 & inlav2obs, ld_nea ) 421 422 CALL obs_mpp_sum_integer( iosdv2obs, iosdv2obsmpp ) 423 CALL obs_mpp_sum_integer( ilanv2obs, ilanv2obsmpp ) 424 CALL obs_mpp_sum_integer( inlav2obs, inlav2obsmpp ) 1065 425 1066 426 ! ----------------------------------------------------------------------- … … 1068 428 ! ----------------------------------------------------------------------- 1069 429 1070 CALL obs_uv_rej( profdata, iuvchku, iuvchkv ) 1071 CALL obs_mpp_sum_integer( iuvchku, iuvchkumpp ) 1072 CALL obs_mpp_sum_integer( iuvchkv, iuvchkvmpp ) 430 IF ( TRIM(profdata%cvars(1)) == 'UVEL' ) THEN 431 CALL obs_uv_rej( profdata, iuvchku, iuvchkv ) 432 CALL obs_mpp_sum_integer( iuvchku, iuvchkumpp ) 433 CALL obs_mpp_sum_integer( iuvchkv, iuvchkvmpp ) 434 ENDIF 1073 435 1074 436 ! ----------------------------------------------------------------------- … … 1109 471 1110 472 IF(lwp) THEN 473 1111 474 WRITE(numout,*) 1112 WRITE(numout,*) 'obs_pre_vel :' 1113 WRITE(numout,*) '~~~~~~~~~~~' 1114 WRITE(numout,*) 1115 WRITE(numout,*) ' Profiles outside time domain = ', & 475 WRITE(numout,*) ' Profiles outside time domain = ', & 1116 476 & iotdobsmpp 1117 WRITE(numout,*) ' Remaining profiles that failed grid search = ', &477 WRITE(numout,*) ' Remaining profiles that failed grid search = ', & 1118 478 & igrdobsmpp 1119 WRITE(numout,*) ' Remaining Udata outside space domain = ', &1120 & iosd uobsmpp1121 WRITE(numout,*) ' Remaining Udata at land points = ', &1122 & ilan uobsmpp479 WRITE(numout,*) ' Remaining '//prodatqc%cvars(1)//' data outside space domain = ', & 480 & iosdv1obsmpp 481 WRITE(numout,*) ' Remaining '//prodatqc%cvars(1)//' data at land points = ', & 482 & ilanv1obsmpp 1123 483 IF (ld_nea) THEN 1124 WRITE(numout,*) ' Remaining Udata near land points (removed) = ',&1125 & inla uobsmpp484 WRITE(numout,*) ' Remaining '//prodatqc%cvars(1)//' data near land points (removed) = ',& 485 & inlav1obsmpp 1126 486 ELSE 1127 WRITE(numout,*) ' Remaining U data near land points (kept) = ',& 1128 & inlauobsmpp 1129 ENDIF 1130 WRITE(numout,*) ' U observation rejected since V rejected = ', & 1131 & iuvchku 1132 WRITE(numout,*) ' U data accepted = ', & 487 WRITE(numout,*) ' Remaining '//prodatqc%cvars(1)//' data near land points (kept) = ',& 488 & inlav1obsmpp 489 ENDIF 490 IF ( TRIM(profdata%cvars(1)) == 'UVEL' ) THEN 491 WRITE(numout,*) ' U observation rejected since V rejected = ', & 492 & iuvchku 493 ENDIF 494 WRITE(numout,*) ' '//prodatqc%cvars(1)//' data accepted = ', & 1133 495 & prodatqc%nvprotmpp(1) 1134 WRITE(numout,*) ' Remaining Vdata outside space domain = ', &1135 & iosdv obsmpp1136 WRITE(numout,*) ' Remaining Vdata at land points = ', &1137 & ilanv obsmpp496 WRITE(numout,*) ' Remaining '//prodatqc%cvars(2)//' data outside space domain = ', & 497 & iosdv2obsmpp 498 WRITE(numout,*) ' Remaining '//prodatqc%cvars(2)//' data at land points = ', & 499 & ilanv2obsmpp 1138 500 IF (ld_nea) THEN 1139 WRITE(numout,*) ' Remaining Vdata near land points (removed) = ',&1140 & inlav obsmpp501 WRITE(numout,*) ' Remaining '//prodatqc%cvars(2)//' data near land points (removed) = ',& 502 & inlav2obsmpp 1141 503 ELSE 1142 WRITE(numout,*) ' Remaining V data near land points (kept) = ',& 1143 & inlavobsmpp 1144 ENDIF 1145 WRITE(numout,*) ' V observation rejected since U rejected = ', & 1146 & iuvchkv 1147 WRITE(numout,*) ' V data accepted = ', & 504 WRITE(numout,*) ' Remaining '//prodatqc%cvars(2)//' data near land points (kept) = ',& 505 & inlav2obsmpp 506 ENDIF 507 IF ( TRIM(profdata%cvars(1)) == 'UVEL' ) THEN 508 WRITE(numout,*) ' V observation rejected since U rejected = ', & 509 & iuvchkv 510 ENDIF 511 WRITE(numout,*) ' '//prodatqc%cvars(2)//' data accepted = ', & 1148 512 & prodatqc%nvprotmpp(2) 1149 513 … … 1151 515 WRITE(numout,*) ' Number of observations per time step :' 1152 516 WRITE(numout,*) 1153 WRITE(numout,997) 517 WRITE(numout,'(10X,A,5X,A,5X,A,A)')'Time step','Profiles', & 518 & ' '//prodatqc%cvars(1)//' ', & 519 & ' '//prodatqc%cvars(2)//' ' 1154 520 WRITE(numout,998) 1155 521 ENDIF … … 1185 551 ENDIF 1186 552 1187 997 FORMAT(10X,'Time step',5X,'Profiles',5X,'Zonal Comp.',5X,'Meridional Comp.')1188 553 998 FORMAT(10X,'---------',5X,'--------',5X,'-----------',5X,'----------------') 1189 554 999 FORMAT(10X,I9,5X,I8,5X,I11,5X,I8) 1190 555 1191 END SUBROUTINE obs_pre_ vel556 END SUBROUTINE obs_pre_prof 1192 557 1193 558 SUBROUTINE obs_coo_tim( kcycle, & … … 1391 756 & kobsno, & 1392 757 & kobsyea, kobsmon, kobsday, kobshou, kobsmin, & 1393 & ktyp, kobsqc, kobsstp, kotdobs, kdailyavtypes, & 1394 & ld_dailyav ) 758 & ktyp, kobsqc, kobsstp, kotdobs, kdailyavtypes ) 1395 759 !!---------------------------------------------------------------------- 1396 760 !! *** ROUTINE obs_coo_tim *** … … 1436 800 INTEGER, DIMENSION(imaxavtypes), OPTIONAL :: & 1437 801 & kdailyavtypes ! Types for daily averages 1438 LOGICAL, OPTIONAL :: ld_dailyav ! All types are daily averages1439 802 !! * Local declarations 1440 803 INTEGER :: jobs … … 1470 833 ENDIF 1471 834 1472 !------------------------------------------------------------------------1473 ! If ld_dailyav is set then all data assumed to be daily averaged1474 !------------------------------------------------------------------------1475 1476 IF ( PRESENT( ld_dailyav) ) THEN1477 IF (ld_dailyav) THEN1478 DO jobs = 1, kobsno1479 1480 IF ( kobsqc(jobs) <= 10 ) THEN1481 1482 IF ( kobsstp(jobs) == (nit000 - 1) ) THEN1483 kobsqc(jobs) = kobsqc(jobs) + 141484 kotdobs = kotdobs + 11485 CYCLE1486 ENDIF1487 1488 ENDIF1489 END DO1490 ENDIF1491 ENDIF1492 835 1493 836 END SUBROUTINE obs_coo_tim_prof … … 1617 960 END DO 1618 961 1619 CALL obs_int_comm_2d( 2, 2, kobsno, igrdi, igrdj, pmask, zgmsk )1620 CALL obs_int_comm_2d( 2, 2, kobsno, igrdi, igrdj, plam, zglam )1621 CALL obs_int_comm_2d( 2, 2, kobsno, igrdi, igrdj, pphi, zgphi )962 CALL obs_int_comm_2d( 2, 2, kobsno, kpi, kpj, igrdi, igrdj, pmask, zgmsk ) 963 CALL obs_int_comm_2d( 2, 2, kobsno, kpi, kpj, igrdi, igrdj, plam, zglam ) 964 CALL obs_int_comm_2d( 2, 2, kobsno, kpi, kpj, igrdi, igrdj, pphi, zgphi ) 1622 965 1623 966 DO jobs = 1, kobsno … … 1713 1056 USE dom_oce, ONLY : & ! Geographical information 1714 1057 & gdepw_1d, & 1715 & gdepw_0, & 1716 #if defined key_vvl 1717 & gdepw_n, & 1718 & gdept_n, & 1719 #endif 1720 & ln_zco, & 1721 & ln_zps 1722 1058 & gdepw_0, & 1059 #if defined key_vvl 1060 & gdepw_n, & 1061 & gdept_n, & 1062 #endif 1063 & ln_zco, & 1064 & ln_zps 1065 1723 1066 !! * Arguments 1724 1067 INTEGER, INTENT(IN) :: kprofno ! Number of profiles … … 1800 1143 END DO 1801 1144 1802 CALL obs_int_comm_3d( 2, 2, kprofno, kp k, igrdi, igrdj, pmask, zgmsk )1803 CALL obs_int_comm_2d( 2, 2, kprofno, igrdi, igrdj, plam, zglam )1804 CALL obs_int_comm_2d( 2, 2, kprofno, igrdi, igrdj, pphi, zgphi )1145 CALL obs_int_comm_3d( 2, 2, kprofno, kpi, kpj, kpk, igrdi, igrdj, pmask, zgmsk ) 1146 CALL obs_int_comm_2d( 2, 2, kprofno, kpi, kpj, igrdi, igrdj, plam, zglam ) 1147 CALL obs_int_comm_2d( 2, 2, kprofno, kpi, kpj, igrdi, igrdj, pphi, zgphi ) 1805 1148 1806 1149 DO jobs = 1, kprofno … … 1827 1170 END DO 1828 1171 1829 ! Check if next to land 1830 IF ( ANY( zgmsk(1:2,1:2,1,jobs) == 0.0_wp ) ) THEN 1831 ll_next_to_land=.TRUE. 1832 ELSE 1833 ll_next_to_land=.FALSE. 1172 ! Check if next to land 1173 IF ( ANY( zgmsk(1:2,1:2,1,jobs) == 0.0_wp ) ) THEN 1174 ll_next_to_land=.TRUE. 1175 ELSE 1176 ll_next_to_land=.FALSE. 1834 1177 ENDIF 1835 1178 1836 1179 ! Reject observations 1837 1180 … … 1850 1193 ENDIF 1851 1194 1852 ! To check if an observations falls within land there are two cases: 1853 ! 1: z-coordibnates, where the check uses the mask 1854 ! 2: terrain following (eg s-coordinates), 1855 ! where we use the depth of the bottom cell to mask observations 1195 ! To check if an observations falls within land there are two cases: 1196 ! 1: z-coordibnates, where the check uses the mask 1197 ! 2: terrain following (eg s-coordinates), 1198 ! where we use the depth of the bottom cell to mask observations 1856 1199 1857 IF( ln_zps .OR. ln_zco ) THEN !(CASE 1) 1858 1859 ! Flag if the observation falls with a model land cell 1860 IF ( SUM( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) & 1861 & == 0.0_wp ) THEN 1862 kobsqc(jobsp) = kobsqc(jobsp) + 12 1863 klanobs = klanobs + 1 1864 CYCLE 1865 ENDIF 1200 IF( ln_zps .OR. ln_zco ) THEN !(CASE 1) 1201 1202 ! Flag if the observation falls with a model land cell 1203 IF ( SUM( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) & 1204 & == 0.0_wp ) THEN 1205 kobsqc(jobsp) = kobsqc(jobsp) + 12 1206 klanobs = klanobs + 1 1207 CYCLE 1208 ENDIF 1866 1209 1867 ! Flag if the observation is close to land 1868 IF ( MINVAL( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) == &1869 & 0.0_wp) THEN 1870 knlaobs = knlaobs + 1 1210 ! Flag if the observation is close to land 1211 IF ( MINVAL( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) == & 1212 & 0.0_wp) THEN 1213 knlaobs = knlaobs + 1 1871 1214 IF (ld_nea) THEN 1872 kobsqc(jobsp) = kobsqc(jobsp) + 14 1873 ENDIF 1874 ENDIF 1215 kobsqc(jobsp) = kobsqc(jobsp) + 14 1216 ENDIF 1217 ENDIF 1875 1218 1876 ELSE ! Case 2 1219 ELSE ! Case 2 1877 1220 1878 ! Flag if the observation is deeper than the bathymetry 1879 ! Or if it is within the mask 1880 IF ( ALL( fsdepw(iig-1:iig+1,ijg-1:ijg+1,kpk) < pobsdep(jobsp) ) & 1881 & .OR. & 1882 & ( SUM( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) & 1883 & == 0.0_wp) ) THEN 1884 kobsqc(jobsp) = kobsqc(jobsp) + 12 1885 klanobs = klanobs + 1 1886 CYCLE 1887 ENDIF 1888 1889 ! Flag if the observation is close to land 1890 IF ( ll_next_to_land ) THEN 1891 knlaobs = knlaobs + 1 1221 ! Flag if the observation is deeper than the bathymetry 1222 ! Or if it is within the mask 1223 IF ( ALL( fsdepw(iig-1:iig+1,ijg-1:ijg+1,kpk) < pobsdep(jobsp) ) & 1224 & .OR. & 1225 & ( SUM( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) & 1226 & == 0.0_wp) ) THEN 1227 kobsqc(jobsp) = kobsqc(jobsp) + 12 1228 klanobs = klanobs + 1 1229 CYCLE 1230 ENDIF 1231 1232 ! Flag if the observation is close to land 1233 IF ( ll_next_to_land ) THEN 1234 knlaobs = knlaobs + 1 1892 1235 IF (ld_nea) THEN 1893 kobsqc(jobsp) = kobsqc(jobsp) + 14 1894 ENDIF 1895 ENDIF 1896 1236 kobsqc(jobsp) = kobsqc(jobsp) + 14 1237 ENDIF 1238 ENDIF 1897 1239 ENDIF 1898 1240 1899 1241 ! For observations on the grid reject them if their are at 1900 1242 ! a masked point … … 1907 1249 ENDIF 1908 1250 ENDIF 1909 1251 1252 ! Flag if the observation falls is close to land 1253 IF ( MINVAL( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) == & 1254 & 0.0_wp) THEN 1255 IF (ld_nea) kobsqc(jobsp) = kobsqc(jobsp) + 14 1256 knlaobs = knlaobs + 1 1257 ENDIF 1910 1258 1911 1259 ! Set observation depth equal to that of the first model depth
Note: See TracChangeset
for help on using the changeset viewer.