GEOS-ESM / GEOSana_GridComp

Repository containing code for data analysis for the GEOS Earth System Model data assimilation
Apache License 2.0
1 stars 3 forks source link

ATMS: this does not look right #172

Closed rtodling closed 1 week ago

rtodling commented 1 week ago

I am looking at QCMOD - for JEDI reasons - and I am finding this logic in qcmod to not make sense for ATMS.

Notice the zeroing out and do loop ranges only make sense if the counter is contiguous. This is incidentally the case for AMSUA, but it is not the case for ATMS.

https://github.com/GEOS-ESM/GEOSana_GridComp/blob/1214e1c0d24eb108ba9861bb8dce6a7a0f08ff28/GEOSaana_GridComp/GSI_GridComp/qcmod.f90#L3118-L3122

The correct thing it seems to me would be to have something like this (the following is cryptic):

ichmw_num = (/ 238, 314, 503, 528, 536, 544, 549, 890 /) if (amuse) ichmw_idx = (/ 1, 2,, 3, 4, 5, 6, 7, 15 /) if (arms) ichmw_idx = (/ 1, 2, 3, 5, 6, 7, 8, 16 /)

then the replacement of the logic:

do i = 1, size(ichmw_idx) errf(i) = zero varinv(i) = zero enddo do i=1,size(ichmw_idx) ii = ichmw_idx(i) if (id_qc(ii) == good_qc) id_qc(ii) * fail_interchan_qc enddo

Other places in QCMOD would have to adapt to the new arrays ichmw_num (not used above but needed elsewhere) and ichmw_idx ...

Needless to say the above takes automate care of the ch 890 handing in the lines that follow the chunk of code above in qcmod.

gmao-yzhu commented 1 week ago

@rtodling The original code should be fine, as additional interchannel QC is done in setuprad.

rtodling commented 1 week ago

@rtodling The original code should be fine, as additional interchannel QC is done in setuprad.

@gmao-yzhu are you saying this is not an error? What I see is that indeed, setupad initializes the quantities being initialized above to zero, but this happens before qc_atms is called; when the gross check is applied, channel index 4 might end up w/ an incorrect errf as far as I can tell. I think there is also something about index 8 too which seem to be handled differently in qc_atms/amsua and setuprad ...

gmao-yzhu commented 1 week ago

@rtodling I think this is not an error. Additional channels of ATMS for interchannel check are handled in setuprad.

gmao-msienkie commented 1 week ago

This code is fine. This code was put in after channel 15 failed in MetOp-B AMSUA, I think. Previously if some lower peaking channels had bad values the entire spot was tossed which would have not allowed MetOpB AMSUA to be assimilated - but this allows the higher peaking channels to be assimilated.

The essence of the test is in the line just above that section you highlighted. If any of the lower peaking channels ( 1,2,3,5,6,7 or 16 for ATMS and 1,2,3,4,5,6 or 16 in AMSUA) has a brightness temperature increment > 200K (probably because of a missing value) then

and then do the reject and fill in ifail_interchan_qc value for the other channels affected for ATMS.

So in essence if any of these key channels (in the 'if' statement above) is missing, reject all the lower peaking channels.

I think that the test was originally written for AMSUA and then "translated" to work for ATMS. The only slightly shady part that I can see is that ATMS channel 4 is not included in the 'if' test. I guess a test of ATMS channel 4 could be added to the 'if' test.

gmao-msienkie commented 1 week ago

OK maybe I should restate this. The channels listed in the if if (any(abs(tbc((/ ich238, ich314, ich503, ich528, ich536, ich544, ich890 /))) &

200.0_r_kind)) then are the channels that are used in subsequent QC tests in the subroutine. So if any of the bc increments for these 'key' channels is wacky we just exclude the low peaking channels - and then all of the various qc tests are in the 'else' section below this 'if'. (And ATMS ch 4 is not relevant in this test.)

EDIT - and when I talk about 'low-peaking' channels I should also clarify that these are the cloud-affected channels that the myriad of tests in this subroutine are supposed to be screening.

rtodling commented 1 week ago

OK maybe I should restate this. The channels listed in the if if (any(abs(tbc((/ ich238, ich314, ich503, ich528, ich536, ich544, ich890 /))) & > 200.0_r_kind)) then are the channels that are used in subsequent QC tests in the subroutine. So if any of the bc increments for these 'key' channels is wacky we just exclude the low peaking channels - and then all of the various qc tests are in the 'else' section below this 'if'. (And ATMS ch 4 is not relevant in this test.)

EDIT - and when I talk about 'low-peaking' channels I should also clarify that these are the cloud-affected channels that the myriad of tests in this subroutine are supposed to be screening.

Hi Meta, thanks for the detailed expiation. That's certainly helpful to understand the context ... however, the explicit list of channels for ATMS - perhaps intentionally - excludes channel 4 and yet, when it comes to zeroing things out in case T > 200 ends up zeroing out the terms related to channel index 4 for ATMS when it seems from the list, this channel should not be included - in any case, the code might be right if its intention is to zero out of info related to the channels w/ index smaller than 8 for ATMS ... but that happens not because of the explicit list of channels indexes but because the do loop is accidentally going through all channels ...

in any case, thank you for the feedback.

gmao-msienkie commented 1 week ago

No, all of the cloud affected channels 1-7, 16, and 16-22 are all rejected and given a zero inverse obs error.
So I can see now that ATMS channel 16 appears to be done twice, once as ich890 and once as a channel in the range 16-22. Channels 17-22 are not participating in the QC tests but they still need to be rejected.

EDIT - you can see there are do loops further down in the code going from 17-22 so probably the loop in this section should be 17-22. But it doesn't really matter that the channel is rejected twice.

rtodling commented 1 week ago

No, all of the cloud affected channels 1-7, 16, and 16-22 are all rejected and given a zero inverse obs error. So I can see now that ATMS channel 16 appears to be done twice, once as ich890 and once as a channel in the range 16-22. Channels 17-22 are not participating in the QC tests but they still need to be rejected.

EDIT - you can see there are do loops further down in the code going from 17-22 so probably the loop in this section should be 17-22. But it doesn't really matter that the channel is rejected twice.

I see all that - I just don't think channel index 4 is handled consistently ... I am closing this .... thanks for the feedback.