- Timestamp:
- 2017-04-18T10:24:44+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/obs_oper_do_not_assim_update/NEMOGCM/NEMO/OPA_SRC/OBS/obs_prep.F90
r7837 r7915 52 52 CONTAINS 53 53 54 SUBROUTINE obs_pre_surf( surfdata, surfdataqc, ld_nea, ld_bound_reject ) 54 SUBROUTINE obs_pre_surf( surfdata, surfdataqc, ld_nea, ld_bound_reject, & 55 kqc_cutoff ) 55 56 !!---------------------------------------------------------------------- 56 57 !! *** ROUTINE obs_pre_sla *** … … 82 83 LOGICAL, INTENT(IN) :: ld_nea ! Switch for rejecting observation near land 83 84 LOGICAL, INTENT(IN) :: ld_bound_reject ! Switch for rejecting obs near the boundary 85 INTEGER, INTENT(IN), OPTIONAL :: kqc_cutoff ! cut off for QC value 84 86 !! * Local declarations 87 INTEGER :: iqc_cutoff = 255 ! cut off for QC value 85 88 INTEGER :: iyea0 ! Initial date 86 89 INTEGER :: imon0 ! - (year, month, day, hour, minute) … … 130 133 ibdysobs = 0 131 134 135 ! Set QC cutoff to optional value if provided 136 IF ( PRESENT(kqc_cutoff) ) iqc_cutoff=kqc_cutoff 137 132 138 ! ----------------------------------------------------------------------- 133 139 ! Find time coordinate for surface data … … 138 144 & surfdata%nsurf, surfdata%nyea, surfdata%nmon, & 139 145 & surfdata%nday, surfdata%nhou, surfdata%nmin, & 140 & surfdata%nqc, surfdata%mstp, iotdobs ) 146 & surfdata%nqc, surfdata%mstp, iotdobs, & 147 & kqc_cutoff = iqc_cutoff ) 141 148 142 149 CALL obs_mpp_sum_integer( iotdobs, iotdobsmpp ) … … 179 186 ALLOCATE( llvalid(surfdata%nsurf) ) 180 187 181 ! We want all data which has qc flags <= 10182 183 llvalid(:) = ( surfdata%nqc(:) <= 10)188 ! We want all data which has qc flags <= iqc_cutoff 189 190 llvalid(:) = ( surfdata%nqc(:) <= iqc_cutoff ) 184 191 185 192 ! The actual copying … … 251 258 & kpi, kpj, kpk, & 252 259 & zmask1, pglam1, pgphi1, zmask2, pglam2, pgphi2, & 253 & ld_nea, ld_bound_reject, kdailyavtypes )260 & ld_nea, ld_bound_reject, kdailyavtypes, kqc_cutoff ) 254 261 255 262 !!---------------------------------------------------------------------- … … 292 299 & pgphi1, & 293 300 & pgphi2 301 INTEGER, INTENT(IN), OPTIONAL :: kqc_cutoff ! cut off for QC value 294 302 295 303 !! * Local declarations 304 INTEGER :: iqc_cutoff = 255 ! cut off for QC value 296 305 INTEGER :: iyea0 ! Initial date 297 306 INTEGER :: imon0 ! - (year, month, day, hour, minute) … … 361 370 iuvchkv = 0 362 371 372 373 ! Set QC cutoff to optional value if provided 374 IF ( PRESENT(kqc_cutoff) ) iqc_cutoff=kqc_cutoff 375 363 376 ! ----------------------------------------------------------------------- 364 377 ! Find time coordinate for profiles … … 371 384 & profdata%nday, profdata%nhou, profdata%nmin, & 372 385 & profdata%ntyp, profdata%nqc, profdata%mstp, & 373 & iotdobs, kdailyavtypes = kdailyavtypes ) 386 & iotdobs, kdailyavtypes = kdailyavtypes, & 387 & kqc_cutoff = iqc_cutoff ) 374 388 ELSE 375 389 CALL obs_coo_tim_prof( icycle, & … … 378 392 & profdata%nday, profdata%nhou, profdata%nmin, & 379 393 & profdata%ntyp, profdata%nqc, profdata%mstp, & 380 & iotdobs )394 & iotdobs, kqc_cutoff = iqc_cutoff ) 381 395 ENDIF 382 396 … … 395 409 396 410 ! ----------------------------------------------------------------------- 397 ! Reject all observations for profiles with nqc > 10398 ! ----------------------------------------------------------------------- 399 400 CALL obs_pro_rej( profdata )411 ! Reject all observations for profiles with nqc > iqc_cutoff 412 ! ----------------------------------------------------------------------- 413 414 CALL obs_pro_rej( profdata, kqc_cutoff = iqc_cutoff ) 401 415 402 416 ! ----------------------------------------------------------------------- … … 419 433 & iosdv1obs, ilanv1obs, & 420 434 & inlav1obs, ld_nea, & 421 & ibdyv1obs, ld_bound_reject ) 435 & ibdyv1obs, ld_bound_reject, & 436 & iqc_cutoff ) 422 437 423 438 CALL obs_mpp_sum_integer( iosdv1obs, iosdv1obsmpp ) … … 440 455 & iosdv2obs, ilanv2obs, & 441 456 & inlav2obs, ld_nea, & 442 & ibdyv2obs, ld_bound_reject ) 457 & ibdyv2obs, ld_bound_reject, & 458 & iqc_cutoff ) 443 459 444 460 CALL obs_mpp_sum_integer( iosdv2obs, iosdv2obsmpp ) … … 452 468 453 469 IF ( TRIM(profdata%cvars(1)) == 'UVEL' ) THEN 454 CALL obs_uv_rej( profdata, iuvchku, iuvchkv )470 CALL obs_uv_rej( profdata, iuvchku, iuvchkv, iqc_cutoff ) 455 471 CALL obs_mpp_sum_integer( iuvchku, iuvchkumpp ) 456 472 CALL obs_mpp_sum_integer( iuvchkv, iuvchkvmpp ) … … 469 485 END DO 470 486 471 ! We want all data which has qc flags = 0472 473 llvalid%luse(:) = ( profdata%nqc(:) <= 10)487 ! We want all data which has qc flags <= iqc_cutoff 488 489 llvalid%luse(:) = ( profdata%nqc(:) <= iqc_cutoff ) 474 490 DO jvar = 1,profdata%nvar 475 llvvalid(jvar)%luse(:) = ( profdata%var(jvar)%nvqc(:) <= 10)491 llvvalid(jvar)%luse(:) = ( profdata%var(jvar)%nvqc(:) <= iqc_cutoff ) 476 492 END DO 477 493 … … 783 799 & kobsno, & 784 800 & kobsyea, kobsmon, kobsday, kobshou, kobsmin, & 785 & ktyp, kobsqc, kobsstp, kotdobs, kdailyavtypes ) 801 & ktyp, kobsqc, kobsstp, kotdobs, kdailyavtypes, & 802 & kqc_cutoff ) 786 803 !!---------------------------------------------------------------------- 787 804 !! *** ROUTINE obs_coo_tim *** … … 827 844 INTEGER, DIMENSION(imaxavtypes), OPTIONAL :: & 828 845 & kdailyavtypes ! Types for daily averages 846 INTEGER, OPTIONAL, INTENT(IN) :: kqc_cutoff ! QC cutoff value 847 829 848 !! * Local declarations 830 849 INTEGER :: jobs 850 INTEGER :: iqc_cutoff=255 831 851 832 852 !----------------------------------------------------------------------- … … 847 867 DO jobs = 1, kobsno 848 868 849 IF ( kobsqc(jobs) <= 10) THEN869 IF ( kobsqc(jobs) <= iqc_cutoff ) THEN 850 870 851 871 IF ( ( kobsstp(jobs) == (nit000 - 1) ).AND.& 852 872 & ( ANY (kdailyavtypes(:) == ktyp(jobs)) ) ) THEN 853 kobsqc(jobs) = kobsqc(jobs) + 14873 kobsqc(jobs) = IBSET(kobsqc(jobs),13) 854 874 kotdobs = kotdobs + 1 855 875 CYCLE … … 894 914 DO jobs = 1, kobsno 895 915 IF ( ( kobsi(jobs) <= 0 ) .AND. ( kobsj(jobs) <= 0 ) ) THEN 896 kobsqc(jobs) = kobsqc(jobs) + 18916 kobsqc(jobs) = IBSET(kobsqc(jobs),12) 897 917 kgrdobs = kgrdobs + 1 898 918 ENDIF … … 906 926 & kobsqc, kosdobs, klanobs, & 907 927 & knlaobs,ld_nea, & 908 & kbdyobs,ld_bound_reject ) 928 & kbdyobs,ld_bound_reject, & 929 & kqc_cutoff ) 909 930 !!---------------------------------------------------------------------- 910 931 !! *** ROUTINE obs_coo_spc_2d *** … … 945 966 LOGICAL, INTENT(IN) :: ld_nea ! Flag observations near land 946 967 LOGICAL, INTENT(IN) :: ld_bound_reject ! Flag observations near open boundary 968 INTEGER, INTENT(IN) :: kqc_cutoff ! Cutoff QC value 969 947 970 !! * Local declarations 948 971 REAL(KIND=wp), DIMENSION(2,2,kobsno) :: & … … 969 992 ! For invalid points use 2,2 970 993 971 IF ( kobsqc(jobs) >= 10) THEN994 IF ( kobsqc(jobs) >= kqc_cutoff ) THEN 972 995 973 996 igrdi(1,1,jobs) = 1 … … 1016 1039 1017 1040 ! Skip bad observations 1018 IF ( kobsqc(jobs) >= 10) CYCLE1041 IF ( kobsqc(jobs) >= kqc_cutoff ) CYCLE 1019 1042 1020 1043 ! Flag if the observation falls outside the model spatial domain … … 1023 1046 & .OR. ( pobsphi(jobs) < -90. ) & 1024 1047 & .OR. ( pobsphi(jobs) > 90. ) ) THEN 1025 kobsqc(jobs) = kobsqc(jobs) + 111048 kobsqc(jobs) = IBSET(kobsqc(jobs),11) 1026 1049 kosdobs = kosdobs + 1 1027 1050 CYCLE … … 1030 1053 ! Flag if the observation falls with a model land cell 1031 1054 IF ( SUM( zgmsk(1:2,1:2,jobs) ) == 0.0_wp ) THEN 1032 kobsqc(jobs) = kobsqc(jobs) + 121055 kobsqc(jobs) = IBSET(kobsqc(jobs),10) 1033 1056 klanobs = klanobs + 1 1034 1057 CYCLE … … 1055 1078 IF (lgridobs) THEN 1056 1079 IF (zgmsk(iig,ijg,jobs) == 0.0_wp ) THEN 1057 kobsqc(jobs) = kobsqc(jobs) + 121080 kobsqc(jobs) = IBSET(kobsqc(jobs),10) 1058 1081 klanobs = klanobs + 1 1059 1082 CYCLE … … 1066 1089 knlaobs = knlaobs + 1 1067 1090 IF (ld_nea) THEN 1068 kobsqc(jobs) = kobsqc(jobs) + 141091 kobsqc(jobs) = IBSET(kobsqc(jobs),9) 1069 1092 CYCLE 1070 1093 ENDIF … … 1101 1124 & kpobsqc, kobsqc, kosdobs, & 1102 1125 & klanobs, knlaobs, ld_nea, & 1103 & kbdyobs, ld_bound_reject ) 1126 & kbdyobs, ld_bound_reject, & 1127 & kqc_cutoff ) 1104 1128 !!---------------------------------------------------------------------- 1105 1129 !! *** ROUTINE obs_coo_spc_3d *** … … 1169 1193 LOGICAL, INTENT(IN) :: ld_nea ! Flag observations near land 1170 1194 LOGICAL, INTENT(IN) :: ld_bound_reject ! Flag observations near open boundary 1195 INTEGER, INTENT(IN) :: kqc_cutoff ! Cutoff QC value 1196 1171 1197 !! * Local declarations 1172 1198 REAL(KIND=wp), DIMENSION(2,2,kpk,kprofno) :: & … … 1196 1222 ! For invalid points use 2,2 1197 1223 1198 IF ( kpobsqc(jobs) >= 10) THEN1224 IF ( kpobsqc(jobs) >= kqc_cutoff ) THEN 1199 1225 1200 1226 igrdi(1,1,jobs) = 1 … … 1246 1272 1247 1273 ! Skip bad profiles 1248 IF ( kpobsqc(jobs) >= 10) CYCLE1274 IF ( kpobsqc(jobs) >= kqc_cutoff ) CYCLE 1249 1275 1250 1276 ! Check if this observation is on a grid point … … 1284 1310 & .OR. ( pobsdep(jobsp) < 0.0 ) & 1285 1311 & .OR. ( pobsdep(jobsp) > gdepw_1d(kpk)) ) THEN 1286 kobsqc(jobsp) = kobsqc(jobsp) + 111312 kobsqc(jobsp) = IBSET(kobsqc(jobsp),11) 1287 1313 kosdobs = kosdobs + 1 1288 1314 CYCLE … … 1299 1325 IF ( SUM( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) & 1300 1326 & == 0.0_wp ) THEN 1301 kobsqc(jobsp) = kobsqc(jobsp) + 121327 kobsqc(jobsp) = IBSET(kobsqc(jobsp),10) 1302 1328 klanobs = klanobs + 1 1303 1329 CYCLE … … 1309 1335 knlaobs = knlaobs + 1 1310 1336 IF (ld_nea) THEN 1311 kobsqc(jobsp) = kobsqc(jobsp) + 141337 kobsqc(jobsp) = IBSET(kobsqc(jobsp),10) 1312 1338 ENDIF 1313 1339 ENDIF … … 1320 1346 & ( SUM( zgmsk(1:2,1:2,kobsk(jobsp)-1:kobsk(jobsp),jobs) ) & 1321 1347 & == 0.0_wp) ) THEN 1322 kobsqc(jobsp) = kobsqc(jobsp) + 121348 kobsqc(jobsp) = IBSET(kobsqc(jobsp),10) 1323 1349 klanobs = klanobs + 1 1324 1350 CYCLE … … 1329 1355 knlaobs = knlaobs + 1 1330 1356 IF (ld_nea) THEN 1331 kobsqc(jobsp) = kobsqc(jobsp) + 141357 kobsqc(jobsp) = IBSET(kobsqc(jobsp),10) 1332 1358 ENDIF 1333 1359 ENDIF … … 1375 1401 END SUBROUTINE obs_coo_spc_3d 1376 1402 1377 SUBROUTINE obs_pro_rej( profdata )1403 SUBROUTINE obs_pro_rej( profdata, kqc_cutoff ) 1378 1404 !!---------------------------------------------------------------------- 1379 1405 !! *** ROUTINE obs_pro_rej *** … … 1393 1419 !! * Arguments 1394 1420 TYPE(obs_prof), INTENT(INOUT) :: profdata ! Profile data 1421 INTEGER, INTENT(IN) :: kqc_cutoff ! QC cutoff value 1422 1395 1423 !! * Local declarations 1396 1424 INTEGER :: jprof … … 1402 1430 DO jprof = 1, profdata%nprof 1403 1431 1404 IF ( profdata%nqc(jprof) > 10) THEN1432 IF ( profdata%nqc(jprof) > kqc_cutoff ) THEN 1405 1433 1406 1434 DO jvar = 1, profdata%nvar … … 1410 1438 1411 1439 profdata%var(jvar)%nvqc(jobs) = & 1412 & profdata%var(jvar)%nvqc(jobs) + 261440 & IBSET(profdata%var(jvar)%nvqc(jobs),14) 1413 1441 1414 1442 END DO … … 1422 1450 END SUBROUTINE obs_pro_rej 1423 1451 1424 SUBROUTINE obs_uv_rej( profdata, knumu, knumv )1452 SUBROUTINE obs_uv_rej( profdata, knumu, knumv, kqc_cutof ) 1425 1453 !!---------------------------------------------------------------------- 1426 1454 !! *** ROUTINE obs_uv_rej *** … … 1442 1470 INTEGER, INTENT(INOUT) :: knumu ! Number of u rejected 1443 1471 INTEGER, INTENT(INOUT) :: knumv ! Number of v rejected 1472 INTEGER, INTENT(IN) :: kqc_cutoff ! QC cutoff value 1473 1444 1474 !! * Local declarations 1445 1475 INTEGER :: jprof … … 1461 1491 DO jobs = profdata%npvsta(jprof,1), profdata%npvend(jprof,1) 1462 1492 1463 IF ( ( profdata%var(1)%nvqc(jobs) > 10) .AND. &1464 & ( profdata%var(2)%nvqc(jobs) <= 10) ) THEN1465 profdata%var(2)%nvqc(jobs) = profdata%var(2)%nvqc(jobs) + 421493 IF ( ( profdata%var(1)%nvqc(jobs) > kqc_cutoff ) .AND. & 1494 & ( profdata%var(2)%nvqc(jobs) <= kqc_cutoff) ) THEN 1495 profdata%var(2)%nvqc(jobs) = IBSET(profdata%var(1)%nvqc(jobs),15) 1466 1496 knumv = knumv + 1 1467 1497 ENDIF 1468 IF ( ( profdata%var(2)%nvqc(jobs) > 10) .AND. &1469 & ( profdata%var(1)%nvqc(jobs) <= 10) ) THEN1470 profdata%var(1)%nvqc(jobs) = profdata%var(1)%nvqc(jobs) + 421498 IF ( ( profdata%var(2)%nvqc(jobs) > kqc_cutoff ) .AND. & 1499 & ( profdata%var(1)%nvqc(jobs) <= kqc_cutoff) ) THEN 1500 profdata%var(1)%nvqc(jobs) = IBSET(profdata%var(1)%nvqc(jobs),15) 1471 1501 knumu = knumu + 1 1472 1502 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
