Changeset 9186 for branches/UKMO/dev_r5518_obs_oper_update_bgc3d/NEMOGCM/NEMO/OPA_SRC/OBS/obs_oper.F90
- Timestamp:
- 2018-01-05T14:29:29+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_obs_oper_update_bgc3d/NEMOGCM/NEMO/OPA_SRC/OBS/obs_oper.F90
r8653 r9186 231 231 DO ji = 1, jpi 232 232 prodatqc%vdmean(ji,jj,jk,1) = 0.0 233 prodatqc%vdmean(ji,jj,jk,2) = 0.0 233 IF ( prodatqc%nvar == 2 ) THEN 234 prodatqc%vdmean(ji,jj,jk,2) = 0.0 235 ENDIF 234 236 END DO 235 237 END DO … … 243 245 prodatqc%vdmean(ji,jj,jk,1) = prodatqc%vdmean(ji,jj,jk,1) & 244 246 & + pvar1(ji,jj,jk) 245 ! Increment field 2 for computing daily mean 246 prodatqc%vdmean(ji,jj,jk,2) = prodatqc%vdmean(ji,jj,jk,2) & 247 & + pvar2(ji,jj,jk) 247 IF ( prodatqc%nvar == 2 ) THEN 248 ! Increment field 2 for computing daily mean 249 prodatqc%vdmean(ji,jj,jk,2) = prodatqc%vdmean(ji,jj,jk,2) & 250 & + pvar2(ji,jj,jk) 251 ENDIF 248 252 END DO 249 253 END DO … … 260 264 prodatqc%vdmean(ji,jj,jk,1) = prodatqc%vdmean(ji,jj,jk,1) & 261 265 & * zdaystp 262 prodatqc%vdmean(ji,jj,jk,2) = prodatqc%vdmean(ji,jj,jk,2) & 263 & * zdaystp 266 IF ( prodatqc%nvar == 2 ) THEN 267 prodatqc%vdmean(ji,jj,jk,2) = prodatqc%vdmean(ji,jj,jk,2) & 268 & * zdaystp 269 ENDIF 264 270 END DO 265 271 END DO … … 297 303 igrdi1(2,2,iobs) = prodatqc%mi(jobs,1) 298 304 igrdj1(2,2,iobs) = prodatqc%mj(jobs,1) 299 igrdi2(1,1,iobs) = prodatqc%mi(jobs,2)-1 300 igrdj2(1,1,iobs) = prodatqc%mj(jobs,2)-1 301 igrdi2(1,2,iobs) = prodatqc%mi(jobs,2)-1 302 igrdj2(1,2,iobs) = prodatqc%mj(jobs,2) 303 igrdi2(2,1,iobs) = prodatqc%mi(jobs,2) 304 igrdj2(2,1,iobs) = prodatqc%mj(jobs,2)-1 305 igrdi2(2,2,iobs) = prodatqc%mi(jobs,2) 306 igrdj2(2,2,iobs) = prodatqc%mj(jobs,2) 305 IF ( prodatqc%nvar == 2 ) THEN 306 igrdi2(1,1,iobs) = prodatqc%mi(jobs,2)-1 307 igrdj2(1,1,iobs) = prodatqc%mj(jobs,2)-1 308 igrdi2(1,2,iobs) = prodatqc%mi(jobs,2)-1 309 igrdj2(1,2,iobs) = prodatqc%mj(jobs,2) 310 igrdi2(2,1,iobs) = prodatqc%mi(jobs,2) 311 igrdj2(2,1,iobs) = prodatqc%mj(jobs,2)-1 312 igrdi2(2,2,iobs) = prodatqc%mi(jobs,2) 313 igrdj2(2,2,iobs) = prodatqc%mj(jobs,2) 314 ENDIF 307 315 END DO 308 316 … … 316 324 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi1, igrdj1, pvar1, zint1 ) 317 325 318 CALL obs_int_comm_2d( 2, 2, ipro, kpi, kpj, igrdi2, igrdj2, plam2, zglam2 ) 319 CALL obs_int_comm_2d( 2, 2, ipro, kpi, kpj, igrdi2, igrdj2, pphi2, zgphi2 ) 320 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi2, igrdj2, pmask2, zmask2 ) 321 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi2, igrdj2, pvar2, zint2 ) 326 IF ( prodatqc%nvar == 2 ) THEN 327 CALL obs_int_comm_2d( 2, 2, ipro, kpi, kpj, igrdi2, igrdj2, plam2, zglam2 ) 328 CALL obs_int_comm_2d( 2, 2, ipro, kpi, kpj, igrdi2, igrdj2, pphi2, zgphi2 ) 329 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi2, igrdj2, pmask2, zmask2 ) 330 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi2, igrdj2, pvar2, zint2 ) 331 ENDIF 322 332 323 333 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi1, igrdj1, pgdept, zgdept ) … … 334 344 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi1, igrdj1, & 335 345 & prodatqc%vdmean(:,:,:,1), zinm1 ) 336 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi2, igrdj2, & 337 & prodatqc%vdmean(:,:,:,2), zinm2 ) 346 IF ( prodatqc%nvar == 2 ) THEN 347 CALL obs_int_comm_3d( 2, 2, ipro, kpi, kpj, kpk, igrdi2, igrdj2, & 348 & prodatqc%vdmean(:,:,:,2), zinm2 ) 349 ENDIF 338 350 339 351 ENDIF … … 378 390 ENDIF 379 391 380 IF ( prodatqc%npvend(jobs,2) > 0 ) THEN 381 382 CALL obs_int_h2d_init( 1, 1, k2dint, zlam, zphi, & 383 & zglam2(:,:,iobs), zgphi2(:,:,iobs), & 384 & zmask2(:,:,1,iobs), zweig2, zmsk_2 ) 385 392 IF ( prodatqc%nvar == 2 ) THEN 393 IF ( prodatqc%npvend(jobs,2) > 0 ) THEN 394 395 CALL obs_int_h2d_init( 1, 1, k2dint, zlam, zphi, & 396 & zglam2(:,:,iobs), zgphi2(:,:,iobs), & 397 & zmask2(:,:,1,iobs), zweig2, zmsk_2 ) 398 399 ENDIF 386 400 ENDIF 387 401 … … 512 526 ENDIF 513 527 514 ! For the second variable 515 IF ( prodatqc%npvend(jobs,2) > 0 ) THEN 516 517 zobsk(:) = obfillflt 518 519 IF ( ANY (idailyavtypes(:) == prodatqc%ntyp(jobs)) ) THEN 520 521 IF ( idayend == 0 ) THEN 522 ! Daily averaged data 528 IF ( prodatqc%nvar == 2 ) THEN 529 ! For the second variable 530 IF ( prodatqc%npvend(jobs,2) > 0 ) THEN 531 532 zobsk(:) = obfillflt 533 534 IF ( ANY (idailyavtypes(:) == prodatqc%ntyp(jobs)) ) THEN 535 536 IF ( idayend == 0 ) THEN 537 ! Daily averaged data 538 539 ! vertically interpolate all 4 corners 540 ista = prodatqc%npvsta(jobs,2) 541 iend = prodatqc%npvend(jobs,2) 542 inum_obs = iend - ista + 1 543 ALLOCATE(interp_corner(2,2,inum_obs),iv_indic(inum_obs)) 544 545 DO iin=1,2 546 DO ijn=1,2 547 548 IF ( k1dint == 1 ) THEN 549 CALL obs_int_z1d_spl( kpk, & 550 & zinm2(iin,ijn,:,iobs), & 551 & zobs2k, zgdept(iin,ijn,:,iobs), & 552 & zmask2(iin,ijn,:,iobs)) 553 ENDIF 554 555 CALL obs_level_search(kpk, & 556 & zgdept(iin,ijn,:,iobs), & 557 & inum_obs, prodatqc%var(2)%vdep(ista:iend), & 558 & iv_indic) 559 560 CALL obs_int_z1d(kpk, iv_indic, k1dint, inum_obs, & 561 & prodatqc%var(2)%vdep(ista:iend), & 562 & zinm2(iin,ijn,:,iobs), & 563 & zobs2k, interp_corner(iin,ijn,:), & 564 & zgdept(iin,ijn,:,iobs), & 565 & zmask2(iin,ijn,:,iobs)) 566 567 ENDDO 568 ENDDO 569 570 ENDIF !idayend 571 572 ELSE 573 574 ! Point data 523 575 524 576 ! vertically interpolate all 4 corners … … 526 578 iend = prodatqc%npvend(jobs,2) 527 579 inum_obs = iend - ista + 1 528 ALLOCATE(interp_corner(2,2,inum_obs),iv_indic(inum_obs)) 529 580 ALLOCATE(interp_corner(2,2,inum_obs), iv_indic(inum_obs)) 581 DO iin=1,2 582 DO ijn=1,2 583 584 IF ( k1dint == 1 ) THEN 585 CALL obs_int_z1d_spl( kpk, & 586 & zint2(iin,ijn,:,iobs),& 587 & zobs2k, zgdept(iin,ijn,:,iobs), & 588 & zmask2(iin,ijn,:,iobs)) 589 590 ENDIF 591 592 CALL obs_level_search(kpk, & 593 & zgdept(iin,ijn,:,iobs),& 594 & inum_obs, prodatqc%var(2)%vdep(ista:iend), & 595 & iv_indic) 596 597 CALL obs_int_z1d(kpk, iv_indic, k1dint, inum_obs, & 598 & prodatqc%var(2)%vdep(ista:iend), & 599 & zint2(iin,ijn,:,iobs), & 600 & zobs2k,interp_corner(iin,ijn,:), & 601 & zgdept(iin,ijn,:,iobs), & 602 & zmask2(iin,ijn,:,iobs) ) 603 604 ENDDO 605 ENDDO 606 607 ENDIF 608 609 !------------------------------------------------------------- 610 ! Compute the horizontal interpolation for every profile level 611 !------------------------------------------------------------- 612 613 DO ikn=1,inum_obs 614 iend=ista+ikn-1 615 616 zweig(:,:,1) = 0._wp 617 618 ! This code forces the horizontal weights to be 619 ! zero IF the observation is below the bottom of the 620 ! corners of the interpolation nodes, Or if it is in 621 ! the mask. This is important for observations near 622 ! steep bathymetry 530 623 DO iin=1,2 531 624 DO ijn=1,2 532 625 533 IF ( k1dint == 1 ) THEN 534 CALL obs_int_z1d_spl( kpk, & 535 & zinm2(iin,ijn,:,iobs), & 536 & zobs2k, zgdept(iin,ijn,:,iobs), & 537 & zmask2(iin,ijn,:,iobs)) 538 ENDIF 539 540 CALL obs_level_search(kpk, & 541 & zgdept(iin,ijn,:,iobs), & 542 & inum_obs, prodatqc%var(2)%vdep(ista:iend), & 543 & iv_indic) 544 545 CALL obs_int_z1d(kpk, iv_indic, k1dint, inum_obs, & 546 & prodatqc%var(2)%vdep(ista:iend), & 547 & zinm2(iin,ijn,:,iobs), & 548 & zobs2k, interp_corner(iin,ijn,:), & 549 & zgdept(iin,ijn,:,iobs), & 550 & zmask2(iin,ijn,:,iobs)) 551 626 depth_loop2: DO ik=kpk,2,-1 627 IF(zmask2(iin,ijn,ik-1,iobs ) > 0.9 )THEN 628 629 zweig(iin,ijn,1) = & 630 & zweig2(iin,ijn,1) * & 631 & MAX( SIGN(1._wp,(zgdepw(iin,ijn,ik,iobs) ) & 632 & - prodatqc%var(2)%vdep(iend)),0._wp) 633 634 EXIT depth_loop2 635 636 ENDIF 637 638 ENDDO depth_loop2 639 552 640 ENDDO 553 641 ENDDO 554 642 555 ENDIF !idayend 556 557 ELSE 558 559 ! Point data 560 561 ! vertically interpolate all 4 corners 562 ista = prodatqc%npvsta(jobs,2) 563 iend = prodatqc%npvend(jobs,2) 564 inum_obs = iend - ista + 1 565 ALLOCATE(interp_corner(2,2,inum_obs), iv_indic(inum_obs)) 566 DO iin=1,2 567 DO ijn=1,2 568 569 IF ( k1dint == 1 ) THEN 570 CALL obs_int_z1d_spl( kpk, & 571 & zint2(iin,ijn,:,iobs),& 572 & zobs2k, zgdept(iin,ijn,:,iobs), & 573 & zmask2(iin,ijn,:,iobs)) 574 575 ENDIF 576 577 CALL obs_level_search(kpk, & 578 & zgdept(iin,ijn,:,iobs),& 579 & inum_obs, prodatqc%var(2)%vdep(ista:iend), & 580 & iv_indic) 581 582 CALL obs_int_z1d(kpk, iv_indic, k1dint, inum_obs, & 583 & prodatqc%var(2)%vdep(ista:iend), & 584 & zint2(iin,ijn,:,iobs), & 585 & zobs2k,interp_corner(iin,ijn,:), & 586 & zgdept(iin,ijn,:,iobs), & 587 & zmask2(iin,ijn,:,iobs) ) 588 589 ENDDO 643 CALL obs_int_h2d( 1, 1, zweig, interp_corner(:,:,ikn), & 644 & prodatqc%var(2)%vmod(iend:iend) ) 645 646 ! Set QC flag for any observations found below the bottom 647 ! needed as the check here is more strict than that in obs_prep 648 IF (sum(zweig) == 0.0_wp) prodatqc%var(2)%nvqc(iend:iend)=4 649 590 650 ENDDO 591 651 652 DEALLOCATE(interp_corner,iv_indic) 653 592 654 ENDIF 593 594 !-------------------------------------------------------------595 ! Compute the horizontal interpolation for every profile level596 !-------------------------------------------------------------597 598 DO ikn=1,inum_obs599 iend=ista+ikn-1600 601 zweig(:,:,1) = 0._wp602 603 ! This code forces the horizontal weights to be604 ! zero IF the observation is below the bottom of the605 ! corners of the interpolation nodes, Or if it is in606 ! the mask. This is important for observations near607 ! steep bathymetry608 DO iin=1,2609 DO ijn=1,2610 611 depth_loop2: DO ik=kpk,2,-1612 IF(zmask2(iin,ijn,ik-1,iobs ) > 0.9 )THEN613 614 zweig(iin,ijn,1) = &615 & zweig2(iin,ijn,1) * &616 & MAX( SIGN(1._wp,(zgdepw(iin,ijn,ik,iobs) ) &617 & - prodatqc%var(2)%vdep(iend)),0._wp)618 619 EXIT depth_loop2620 621 ENDIF622 623 ENDDO depth_loop2624 625 ENDDO626 ENDDO627 628 CALL obs_int_h2d( 1, 1, zweig, interp_corner(:,:,ikn), &629 & prodatqc%var(2)%vmod(iend:iend) )630 631 ! Set QC flag for any observations found below the bottom632 ! needed as the check here is more strict than that in obs_prep633 IF (sum(zweig) == 0.0_wp) prodatqc%var(2)%nvqc(iend:iend)=4634 635 ENDDO636 637 DEALLOCATE(interp_corner,iv_indic)638 639 655 ENDIF 640 656
Note: See TracChangeset
for help on using the changeset viewer.