Skip to content

Feature/update for viirs da #393

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
rtodling merged 5 commits into from
Jan 14, 2026
Merged

Feature/update for viirs da #393

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
c5dd020
starting to add VIIRS data
vbuchard Aug 20, 2025
1695176
moving to openMP QC with no issue
vbuchard Sep 16, 2025
571371f
update in the name of file for NOAA20
vbuchard Dec 3, 2025
b4dbdee
remove SNPP filename for now
vbuchard Dec 3, 2025
b9480eb
Merge branch 'develop' into feature/update_for_VIIRS_DA
rtodling Jan 14, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions src/Applications/GAAS_App/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,7 @@ esma_add_library(${this}
SRCS ${SRCS}
DEPENDENCIES GMAO_psas GMAO_ods GMAO_hermes Chem_Base)

target_link_libraries(${this} PRIVATE OpenMP::OpenMP_Fortran)

ecbuild_add_executable(TARGET ana_aod.x SOURCES ana_aod.F LIBS ${this} GMAO_gfio_r4 Chem_Base Chem_Shared)
ecbuild_add_executable(TARGET mpana_aod.x SOURCES mpana_aod.F90 LIBS ${this} GMAO_gfio_r4 Chem_Base Chem_Shared)
Expand Down
7 changes: 7 additions & 0 deletions src/Applications/GAAS_App/ana.rc.tmpl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -37,10 +37,17 @@ ODS_files::
#___AQUA___${FVWORK}/nnr_003.MYD04_L2a.ocean.%y4%m2%d2_%h200z.ods
#___AQUA___${FVWORK}/nnr_003.MYD04_L2a.deep.%y4%m2%d2_%h200z.ods

#___NOAA20___${FVWORK}/nnr_001.NOAA20_L2a.dt_land.%y4%m2%d2_%h200z.ods
#___NOAA20___${FVWORK}/nnr_001.NOAA20_L2a.dt_ocean.%y4%m2%d2_%h200z.ods
#___NOAA20___${FVWORK}/nnr_001.NOAA20_L2a.db_deep.%y4%m2%d2_%h200z.ods


#___MISR___${FVWORK}/misr_F12_0022.bright_tc8.obs.%y4%m2%d2.ods

#___AERONET___${FVWORK}/aeronet.obs.%y4%m2%d2_%h2z.ods

#___LUNAR_AERONET___${FVWORK}/aeronet.obs.%y4%m2%d2_%h2z.ods

# Passive data
#/nobackup/3/PARASOL/Level2/ODS/Y%y4/M%m2/PARASOL_L2.aero_tc8.obs.%y4%m2%d2.ods
#/nobackup/3/OMI/Level2/ODS/Y%y4/M%m2/omi.aero_tc8.obs.%y4%m2%d2.ods
Expand Down
75 changes: 50 additions & 25 deletions src/Applications/GAAS_App/m_sqc.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -148,7 +148,7 @@ subroutine SQC ( y, ntotal, npassed )
!
! 15May2002 Dee Created this module from psas_qc.f version 1.61
! 31Mar2003 (Dee/Rukh) Minor changes to the buddy check
! 13Oct2020 (Todling) Store sigO in xvec
! 13Oct2020 (Todling) Store sigO in xvec!
!
!EOP
!-------------------------------------------------------------------------
Expand Down Expand Up @@ -801,7 +801,7 @@ subroutine Get_ErrVar_( nobs, nexcl, y, varF, varO )
ktPSAS(np) = kta
end if
end do

! Move those data up...
! ---------------------
nmoved = 0
Expand Down Expand Up @@ -1267,6 +1267,7 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )
! Index vector for sorts (cannot be allocatable because of OMP)
! -------------------------------------------------------------
integer indx(nobs)
integer, allocatable :: indx2(:)

! List of suspect observations
! ----------------------------
Expand All @@ -1275,8 +1276,9 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )

! List of buddies
! ---------------
integer :: ki_buddy(nobs) ! index of buddy
real*8 :: wt_buddy(nobs) ! weight of buddy
integer, allocatable :: ki_buddy(:) ! index of buddy
real*8, allocatable :: wt_buddy(:) ! weight of buddy
real, allocatable :: mysep(:,:)

! Suspect and reaccept flags
! --------------------------
Expand Down Expand Up @@ -1304,7 +1306,7 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )
integer rc, iter, ireg, ibeg, iend, ilen, i, is, ib, ibb, j
integer kis, kts, krs, niter
integer maxreg ! number of PSAS regions
real*8 tol2, exponent, scgain
real*8 tol2, myexponent, scgain

real*8 dist2, z_dist2

Expand Down Expand Up @@ -1448,16 +1450,37 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )
' with ', n_susp, ' suspect observations'

n_reacc = 0
allocate(mysep(n_susp, maxreg))

!!!$omp parallel do default(none) private(is, krs) shared(n_susp, maxreg, kr_susp, partition, mysep)
do is = 1, n_susp
krs = kr_susp(is)
do ireg = 1, maxreg
mysep(is,ireg) = Separation(GetRegion(ireg,partition),GetRegion(krs,partition),SEPANG_MIN)
end do
end do
!!!$omp end parallel do


!$omp parallel do &
!$omp default(shared), &
!$omp private(is,kis,kts,krs,nbuddy,ireg,ibeg,iend,i,exponent, &
!$omp scgain,ki_buddy,wt_buddy,indx,accum_del,accum_de2, &
!$omp accum_wgt,accum_var,ib,ibb,del_star,alpha,tol2), &
!$omp reduction(+:n_reacc)
!$omp& schedule(dynamic) &
!$omp& default(none) &
!$omp& firstprivate(nobs) &
!$omp& private(lvs,is,kis,kts,krs,nbuddy,ireg,ibeg,iend,i,myexponent, &
!$omp& scgain,accum_del,accum_de2, ki_buddy, wt_buddy, indx2,&
!$omp& accum_wgt,accum_var,ib,ibb,del_star,alpha,tol2) &
!$omp& shared(n_susp, ki_susp, kt, kr_susp, lev, maxreg, &
!$omp& ireglen, seplim, iregbeg, single_level, issuspect, &
!$omp& qcx, xobs, yobs, zobs, ls_h, ls_v, search_rad, &
!$omp& varF, VarO, nbuddy_max, OmF, nstar, tau_buddy, reaccept,mysep) &
!$omp& reduction(+:n_reacc)

do is = 1, n_susp ! for each suspect obs

do is = 1, n_susp ! for each suspect obs
allocate(ki_buddy(nobs))
allocate(wt_buddy(nobs))
allocate(indx2(nobs))

kis = ki_susp(is) ! this suspect's index
kts = kt(kis) ! this suspect's kt
krs = kr_susp(is) ! this suspect's region
Expand All @@ -1468,9 +1491,8 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )
nbuddy = 0
do ireg = 1, maxreg

if ((Separation(GetRegion(ireg,partition),GetRegion(krs,partition),SEPANG_MIN) <= seplim) &
if ( (mysep(is,ireg) <= seplim) &
.and. (ireglen(ireg) >0 ) ) then ! nearby region with data

ibeg = iregbeg(ireg)
iend = ibeg + ireglen(ireg) - 1

Expand All @@ -1484,15 +1506,15 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )
if ( .not. issuspect(i) .AND. qcx(i)==0 .AND. &
kt(i)==kts) then

exponent = dist2(kis,i)/ls_h(kts)**2
myexponent = dist2(kis,i)/ls_h(kts)**2

if ( ls_v(kts)>tol_rel ) then ! upper-air data
exponent = exponent + z_dist2(kis,i)/ls_v(kts)**2
myexponent = myexponent + z_dist2(kis,i)/ls_v(kts)**2
end if

! only if not too far (search_rad length scales):

if ( exponent<search_rad**2 ) then
if ( myexponent<search_rad**2 ) then

! found a candidate buddy:

Expand All @@ -1502,7 +1524,7 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )
! associate weight with this candidate:

scgain = varF(i) / (varF(i) + varO(i))
wt_buddy(nbuddy) = scgain * exp(-0.5*exponent)
wt_buddy(nbuddy) = scgain * exp(-0.5*myexponent)

end if ! not too far

Expand All @@ -1518,9 +1540,9 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )

! Find the best buddies
! ---------------------
call IndexSet ( nbuddy, indx )
call IndexSort ( nbuddy, indx, wt_buddy, descend=.true. )


call IndexSet ( nbuddy, indx2 )
call IndexSort ( nbuddy, indx2, wt_buddy, descend=.true. )
nbuddy = min(nbuddy, nbuddy_max) ! pick highest weights

! Do buddy check
Expand All @@ -1532,14 +1554,13 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )

do ib = 1, nbuddy ! accumulate weights, variances, etc.

ibb = indx(ib) ! index in ki_buddy of best buddy
ibb = indx2(ib) ! index in ki_buddy of best buddy
i = ki_buddy(ibb) ! index in OmF of best buddy

accum_del = accum_del + OmF(i) * wt_buddy(ibb)
accum_wgt = accum_wgt + wt_buddy(ibb)
accum_de2 = accum_de2 + (OmF(i))**2
accum_var = accum_var + varO(i) + varF(i)

end do

del_star = sngl(accum_del/accum_wgt) ! prediction
Expand All @@ -1560,10 +1581,15 @@ subroutine Buddy_Check_ ( nobs, nfailed, y, varF, varO )

end if

deallocate(indx2)
deallocate(ki_buddy)
deallocate(wt_buddy)

end do

!$omp end parallel do

deallocate(mysep)
if (verbose) then
write(stdout,'(2a,i8,a,i2)') myname, ': ', &
n_reacc,' observations reaccepted after iteration', iter
Expand Down Expand Up @@ -1629,7 +1655,6 @@ subroutine Region_Sort_ ( )
write(stderr,'(2a)') myname,': allocate() error'
call die(myname)
end if

! Compute iregn from (x,y,z)
! -------------------------
call XYZ2REG ( nobs, xobs, yobs, zobs, kr, partition )
Expand All @@ -1644,9 +1669,9 @@ subroutine Region_Sort_ ( )
zobs = zobs( (/ (indx(i), i=1,nobs) /) )
varF = varF( (/ (indx(i), i=1,nobs) /) )
varO = varO( (/ (indx(i), i=1,nobs) /) )

! Also re-sort local arrays
! -------------------------
! --------------------------
lon = lon( (/ (indx(i), i=1,nobs) /) )
lat = lat( (/ (indx(i), i=1,nobs) /) )
OmF = OmF( (/ (indx(i), i=1,nobs) /) )
Expand Down
52 changes: 52 additions & 0 deletions src/Applications/GAAS_App/psas.rc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -129,6 +129,18 @@ DataSourceTable::
320 MYDD -100 MODIS/AQUA-Land (Deep Blue)
323 ANET -100 AERONET
324 AVHRR -100 AVHRR NNR Retrievals
335 VSNPPDBD -100 VIIRS SNPP Deep Blue Deep
336 VSNPPDTL -100 VIIRS SNPP Dark Target Land
337 VSNPPDTO -100 VIIRS SNPP DARK Target Ocean
338 VN20DBL -100 VIIRS NOAA20 Deep Blue LAND
339 VN20DBO -100 VIIRS NOAA20 Deep Blue OCEAN
340 VN20DBD -100 VIIRS NOAA20 Deep Blue DEEP
341 VN20DTL -100 VIIRS NOAA20 Dark Target LAND
342 VN20DTO -100 VIIRS NOAA20 Dark Target OCEAN
343 LANET -100 Lunar Aeronet
344 VN21DTL -100 VIIRS NOAA21 Dark Target LAND
345 VN21DTO -100 VIIRS NOAA21 Dark Target OCEAN
346 VN21DBD -100 VIIRS NOAA21 Deep Blue DEEP

:: # End of DataSourceTable

Expand Down Expand Up @@ -206,6 +218,46 @@ ObsErr*AVHRR:: # kx = 324 (like MYDO)
q_UprAir.u 0.22 0.21 0.20 0.19
::

ObsErr*VSNPPDBD:: # kx = 335
q_UprAir.u 0.22 0.21 0.20 0.19
::

ObsErr*VSNPPDTO:: # kx = 336
q_UprAir.u 0.19 0.19 0.18 0.18
::

ObsErr*VSNPPDTL:: # kx = 337
q_UprAir.u 0.19 0.19 0.18 0.18
::

ObsErr*VN20DTL:: # kx = 341
q_UprAir.u 0.19 0.19 0.18 0.18
::

ObsErr*VN20DTO:: # kx = 342
q_UprAir.u 0.19 0.19 0.18 0.18
::

ObsErr*VN20DBD:: # kx = 340
q_UprAir.u 0.22 0.21 0.20 0.19
::

ObsErr*VN21DTL:: # kx = 344
q_UprAir.u 0.19 0.19 0.18 0.18
::

ObsErr*VN21DTO:: # kx = 345
q_UprAir.u 0.19 0.19 0.18 0.18
::

ObsErr*VN21DBD:: # kx = 346
q_UprAir.u 0.22 0.21 0.20 0.19
::

ObsErr*LANET:: # kx = 343 -- same as ANET for now
q_UprAir.u 0.22 0.21 0.20 0.19
::

#...............................................................................

# --------------------------------------------------
Expand Down