1 / 12

CrIMSS EDR offline QC Report Response to recent findings at NOAA STAR

CrIMSS EDR offline QC Report Response to recent findings at NOAA STAR. Susan Kizer and Xu Liu NASA Langley Research Center Murty Divakarla , Mike Wilson, Xiaozhen Xiong , Changyi Tan, and Flavio Iturbide IM Systems Group, Inc. at NOAA/STAR

Télécharger la présentation

CrIMSS EDR offline QC Report Response to recent findings at NOAA STAR

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript


  1. CrIMSS EDR offline QC Report Response to recent findings at NOAA STAR Susan Kizer andXu Liu NASA Langley Research Center MurtyDivakarla, Mike Wilson, XiaozhenXiong, ChangyiTan, and Flavio IturbideIM Systems Group, Inc. at NOAA/STAR DeguiGu, Denise Hagan, and Xia L Ma Northrop Grumman Aerospace Systems Eric Maddy, Antonia Gambacorta, Nick Nalli, Bigyani Das IM Systems Group, Inc. at NOAA/STAR Chris Barnet&, Tony Reale*, and Mitch Goldberg$&Formerly with NOAA/STAR and currently with STC, Columbia, MD*NOAA/STAR $JPSS/NOAA

  2. Today’s Discussion • Recent Findings from Mike Wilson and Bigyani Das regarding a bug in assigning quality control flags • LaRCproposed solutions to these findings • Discrepancy Report submitted by Xu Liu • Initialization of ccNAF • Relevant to this issue • Description of Quality Control flags reported to netCDF output file • Use IR+MW flags for an example • Mention new CDFCB released Nov/Dec 2012 with updated overall quality QC flag

  3. NOAA STAR QC Discrepancy Findings Bigyani found an issue in my qc chart code that showed different yields when you look at the QC flags, versus when you generate them by hand (using chi-square boundaries).  I tracked down the discrepancy, and it turns out to be a bug in the EDR code.  Here's a quick sample from my program:Category Perc.  QC(1) QC(4)  QC(5) QC(13)---------------------------------------------------------------------------------------------------------------------------------------------------All Categories                                            100.00 77.07 32.1886.3330.81Combined Retrieval Converged                          31.62 97.44 99.31 96.85 97.44Combined Retrieval Did Not Converge              68.38 67.66 1.13 81.47 0.00MW-Only Retrieval Converged                            86.70 82.59 36.18 99.58 34.70MW-Only Retrieval Did Not Converge            13.30 41.15 6.11 0.00 5.47The numbers in red should be an exact match to each other.  The numbers in blue should also match each other.  They don't.  Here's the difference:For the MW-only run:86.33%: # of profiles with a chi-square value less than or equal to 2.86.70%: # of profiles with a passing QC flag for the MW-only run.For the combined run:30.81%: # of profiles with a MW chi-square <= 2 (looking at MX6.6 here) and IR chi-square <= 131.62%: # of profiles with a passing QC flag for the MW+IR run.These should be the same value, just calculated different ways. 

  4. Here's what's happening.In the calcCrimssProfiles.f code, you initialize three arrays with 999999 to hold chisq values:chisq_ir = CHISQ_INIT_VALUE             chisq_mw1 = CHISQ_INIT_VALUE             chisq_mw2 = CHISQ_INIT_VALUE We don't actually use these when we're calculated chisq values.  Instead we have "working" values called chisq, chisqmw1, and chisqmw.  These values are inside the ir loop and have the chisq for the current iteration.  At the end of each ir+mw iteration, you write these values to chisq_ir(iter), chisq_mw1(iter), and chisq_mw2(iter) respectively.  When the loop is finished and you have a final answer, you call qc using chisq and chisqmw2.  That gets reported.In a few special cases, this falls apart.  When you have overcast skies, you have this code in calcCrimssProfiles.f:               ELSE !                 xma (11/15/06) skip below part and write out MW retrievals                   IF (writeDebug) THEN CALL addLine('Overcast case, write out MW retrievals')                   END IF qcflag(1) = 1 qcflag(4) = 1                   EXIT irLoop It looks good.  But the call to qc is after the loop ends:CALL qc(chisq,chisqmw,chisqairs,nsurf,xgesg,xgesmwsav,plandavg, &qcflag, iflagqc) We're not passing our initialized variables in.  Instead, we're sending in chisq and chisqmw which are likely left over from the last non-overcast profile.  Our values for qcflag(1) and qcflag(4) get overwritten by whatever happens in qc, and now the flags we set to 1 are often reset to 0.There's another case where this happens.  If ccNaf went too high, we have this bit in the code:               IF (( .NOT. clrflg(ic)) .AND. iter == maxIrIter .AND. & sqrt(ccnaf)>5.) THEN chisq_ir(iter) = 100000.0 + chisq_ir(iter)                   EXIT irLoop                END IF Again, the call to qc is using chisq instead of chisq_ir(iter).  So the qc routine never sees the 100000 added on, and some values pass. In the case of the microwave-only product QC(5), I believe it's the same issue causing the reverse problem.  The mw-only QC pass/fail is set early, and that value seems right.  What happens later is this:            chisq_mw1 = CHISQ_INIT_VALUEThe final chisq value that determines QC(5) is stored in chisqmw1.  It then gets reassigned here:               chisq_mw1(iter) = chisqmw1 (This is also why we report the same chisq value up to 4 times when we look up the value of mw1 in the final product.  Different issue.)  In cases where we add the 100000, we're fine.   But when we have overcast skies, the line for chisq_mw1 is never reached.  So we leave chisq_mw1 as a 999999 and that gets reported.  But the QC may either pass or fail.One possible fix is to do this before qc?chisq = chisq_ir(iter)chisqmw = chisq_mw2(iter)                               chisq_mw1(iter)=chisqmw1

  5. LaRC Response to QC Discrepancy • In the most basic terms, the discrepancies are introduced in the instances of “EXIT irLoop”. • LaRC proposed solutions: • Replace the initialization of chisq_mw1 to FILL by setting it equal to the 1st stage MW chisq value. Remove this assignment from within the IR+MW iteration. • chisq_mw1(1:maxIrIter) = chisqmw1 irLoop: DO iter = 1,mxIrIter • Do not use the “working” variables for chisq in the QC subroutine. Instead pass chisq_ir and chisq_mw2 CALL qc(chisq_ir(iter),chisq_mw2(iter), chisqairs,nsurf,xgesg,xgesmwsav,plandavg,qcflag,iflagqc) • Also, after the 1st stage MW retrieval, chisqmw is stored in the structure that is passed to the post-processor (to create the EDR). Currently, that parameter is set to a variable that has not yet been assigned crimssRetrLvlDataPtr%chisqMw(currRet) = chisq_mw1(iter) I propose that this line be removed? The same line occurs after the IR+MW retrieval after chisq_mw1 has been assigned.

  6. ccNAF Discrepancy Report Issue: Currently, the IDPS code does not give ccNAF proper values under 2 conditions: 1). When the scene is identified as clear 2). When the algorithm skipped combined IR+MW retrieval Therefore, as it stands now, ccNAF is not a good QC parameter for the current IDPS code. Also, it is important to note that ccNAF is not the true noise amplification factor. It is the square of the NAF. It is passed to the c++ data structure by the following line in “calcCrimssProfiles.f”. crimssRetrLvlDataPtr%IrNAFactor(currRet) = sqrt(ccnaf) It can be easily written to the IP product file if it is not already written. For the off-line code, we do write it out. For the idps code, it passes back to DMS and we have not trace back that far.

  7. ccNAF Discrepancy Report Proposed solution: 1). Initialize ccNAF to 999999.0 (or something appropriate). Add ccNAF= 999999.0 before these two lines in the IDPS code “calcCrimssProfiles.f” ! Start IR iteration: irLoop: DO iter = 1, maxIrIter 2). Assign ccNAF if doing clear-sky retrieval. Add ccNAF=1.0/9 after the following code segment in “calcCrimssProfiles.f”: IF (clrflg(ic)) THEN IF (iter <= 1) THEN CALL avgrad(ym0,ic,nfovcc, & indxfov,kchan,kchanFOV,devnoise,ym,ef) 3). Comment out the following line in the “calc_cc_rad.f” IF (ccnaf<1.0) ccnaf = 1.0 (This change may affect the CC results. We may need to do a performance testing before committing this change)

  8. Description of netCDF QC flagsExample taken from IR+MW file qcFlagIr(16) = 1 MW Only Retrieval - bug identified in reporting this flag, unreliable currently qcFlagIr(17) = 1 no convergence - bug identified in reporting this flag, unreliable currently qcFlagIr(18) = 0 IR-MW profile diff within threshold, = 1 IR-MW profile diff exceeded threshold qcFlagIr(19) = 0 Clear Scene, = 1 Partly Cloudy Scene qcFlagIr(20) = 0 Local Thermodynamic Equilibrium, = 1 Non- Local Thermodynamic Equilibrium qcFlagIr(21) = 0 No precipitation, = 1 Precipitation Detected qcFlagIr(22) = 0 9 FOVs used, = 1 4 FOVs used - bug identified in reporting this flag, unreliable currently qcFlagIr(23) = 0 9 FOVs used, = 1 1 FOV used - bug identified in reporting this flag, unreliable currently qcFlagIr(24) = 0 9 FOVs used, = 1 No retrieval FOV - bug identified in reporting this flag, unreliable currently qcFlagIr(25) = 0 Temperature within range, = 1 Temperature out of range qcFlagIr(26) = 0 no sun glint, = 1 sun glint qcFlagIr(27) = 1 ATMS not available qcFlagIr(28) = 0 daytime, = 1 nighttime qcFlagIr(29) = 0 no ice, = 1 ice mask qcFlagIr(30) = ascending/descending flag - currently not set qcFlagIr(31) = 1 ATMS quality high Quality Control Flags(nProfiles, 31) qcFlagIr(1) = 1 IR+MW retrieval chisq not within threshold qcFlagIr(2) = 1, IR-MW profile difference not within threshold (over ocean) qcFlagIr(3) = 1 chisqAirs not within threshold qcFlagIr(4) = 1 2nd stage MW chisq not within threshold qcFlagIr(5) = 1 MW did not converge qcFlagIr(6) = 1 precipitation detected qcFlagIr(7) = 1 MW not available qcFlagIr(8) = 0 day, = 1 night qcFlagIr(9) = 0 clear scene, = 1 partly cloudy scene, = 2 cloudy scene qcFlagIr(10) = 0 ocean, = 1 land, = 2 coast qcFlagIr(11) = 1 ice over ocean qcFlagIr(12) = 1 non-LTE condition qcFlagIr(13) = 0 IR+MW retrieval converged, = 1 IR+MW retrieval did not converge qcFlagIr(14) = 0 MW-only retrieval converged, = 1 MW-only retrieval did not converge qcFlagIr(15) = 1 IR+MW Only Retrieval - bug identified in reporting this flag, unreliable currently

  9. NOTES regarding QC Flag IssueMW retrieval IF MW data available THEN Stage 1 MW Retrieval • Setup radiances and noise covariances • 100 CONTINUE • Setup initial guess profile and Generate radiance based upon initial guess, calculate residual • Calculate closeness between simulated radiance and initial first guess, chisqmw0 • Start MW iteration chisq_mw = -999.5 Added in offline code, initialization for netCDF output. • DO itermw = 1, maxMwIter • Generate MW radiance • Check for convergence, chisqmw • chisq_mw(itermw) = chisqmwNot in IDPS code, for netCDF output. • Stratification scheme: reassign cov based on MW retrieval • END DO • chisqmw1 = chisqmw • Set the MW only retrieval not converged bit flag if QC threshold was exceeded qcflag(5) = 1 ELSE MW Data not available END IF crimssRetrLvlDataPtr%chisqMw(currRet) = chisq_mw1(iter) chisq_mw1 has not yet been assigned a value . . . Output MW to netCDF End of stage 1 MW Retrieval

  10. NOTES regarding QC Flag IssueIR+MW Retrieval Proceed with IR retrieval Generate MW and IR radiance . . . Initialize variables used during iteration chisq0 = CHISQ_INIT_VALUE chisq_ir = CHISQ_INIT_VALUE chisq_mw1 = CHISQ_INIT_VALUE chisq_mw2 = CHISQ_INIT_VALUE Start IR iteration: irLoop: DO iter = 1, maxIrIter Propose to move here IF (cc(ic)) THEN chisq_mw1(iter) = chisqmw1 IF (clrflg(ic)) THEN IF (iter <= 1) THEN . . . If all channels are turned off EXIT irLoopchisq has not been initialized END IF (used in defining QC flags in CALL to qc) NOTE: Assuming this is independent of iteration and will occur when iter = 1. Therefore chisq_ir, etc., will still be initialized to CHISQ_INIT_VALUE when passed to CALL qc. chisq_mw1 not yet assigned a value so it will be reported as FILL (999999).

  11. NOTES regarding QC Flag IssueIR+MW Retrieval ELSE Generate cloud cleared radiance, assign value to ccNAF cc(ic) = TRUE, clrflg(ic) = FALSE END IF ELSE Overcast Case, set qcflag(1) = 1 and qcflag(4) = 1 EXIT irLoopchisq has not been initialized; END IF qcflag(1) and qcflag(4) will be reset to = 0 NOTE: Assuming this is independent of iteration and will occur when iter = 1. Therefore chisq_ir, etc., will still be initialized to CHISQ_INIT_VALUE (999999) when passed to CALL qc. chisq_mw1 not yet assigned a value so it will be reported as FILL (999999). Generate updated radiances based upon latest retrieval Check convergence criteria: chisqmw, chisq, qcchisq chisq_ir(iter) = chisq chisq_mw1(iter) = chisqmw1 Propose to move this line up to beginning of chisq_mw2(iter) = chisqmwirLoop. IF (.NOT. cc(ic)) THEN This will not be reached. ccnaf = 0.0001 END IF

  12. NOTES regarding QC Flag IssueIR+MW Retrieval … IF ((.NOT. Clrflg(ic)) AND iter == maxIrIter AND sqrt(ccNAF) > 5.) THEN ccNAF not initialized. Value set only when chisq_ir(iter) = 100000.0 + chisq_ir(iter) cc(ic) = TRUE AND clrflg(ic) = FALSE EXIT irLoop Needs to have impact on CALL to qc END IF but will not if chisq_ir is not passed in. … END DO irLoop Set the chi-square values for the current retrieval crimssRetrLvlDataPtr%chisq(currRet) = chisq_ir(iter) crimssRetrLvlDataPtr%chisqMw(currRet) = chisq_mw1(iter) crimssRetrLvlDataPtr%chisqMw2(currRet) = chisq_mw2(iter) … crimssRetrLvlDataPtr%IrNAFactor(currRet) = sqrt(ccnaf) ... CALL qc(chisq, chisqmw, chisqairs, nsurf, xgesg, xgesmwsav, plandavg, Propose changing to qcflag, iflagqc) CALL qc(chisq_ir(iter), chisq_mw2(iter), … This subroutine assigns values to qcflag(1:4) based on chisq values, profile differences between MW1 and MW2, and land fraction. … IF IR+MW retrieval good, save it as 1st guess for the next cluster IF (chisq < 2. AND chisqmw < 2.) THEN Is it best to make this chisq_ir(iter) and … chisq_mw2(iter) as well? END IF

More Related