channel subsetting resulting in potentially different TB values

[BenjaminTJohnson]

Reported by Yaping Zhou, 8/13/2024 (possibly only applies to v2.3, will verify in v2.x and v3): (edit: reported as v2.1.1)

"
I have been trying to run a subset of channels in CRTM with MiRS channel selection since this will save forward model running time in the retrieval.

The version of CRTM within the MiRS system has already implemented mechanism for running only selected channels by using the ChannelInfo%Process_Channel flag but MiRS does not have mechanism to invoke that function. So it’s a matter I pass the channel selection information to that Process_Channel parameter. I am mostly successful in passing that information to CRTM and only run a subset in the channel_loop in CRTM_K_Matrix.f90. However, I noticed that the TB returned from the same channel is different if the channel index ln is accumulated using effectively computed channels or accumulated with all channels:

 Channel_Loop: DO l = 1, ChannelInfo(n)%n_Channels
          ln = ln + 1                                   !yp include all channels whether calculated or not

          ! Channel setup
          ! ...Skip channel if requested
          !IF ( .NOT. ChannelInfo(n)%Process_Channel(l) ) print*, 'skip channel ', l 
          IF ( .NOT. ChannelInfo(n)%Process_Channel(l) ) CYCLE Channel_Loop
 
          ! ...Shorter name
          ChannelIndex = ChannelInfo(n)%Channel_Index(l)

          ! ...Increment the processed channel counter
          !ln = ln + 1                                        !yp only include computed channels

          ! ...Assign sensor+channel information to output 

After TB calculation, I printed out the BT for the first channel I selected (no. 12 channel) in the second sensor. My first sensor is ATMS, I used all of the channels. So when ln is accumulated after cycle channel_loop, ln for this channel is 23, while if ln is accumulated before cycle channel_loop, ln is 22+12=34, and here is the printout from the two situations:

           CALL CRTM_Planck_Temperature( &
                   SensorIndex                            , & ! Input
                   ChannelIndex                           , & ! Input
                   RTSolution(ln,m)%Radiance              , & ! Input
                   RTSolution(ln,m)%Brightness_Temperature  ) ! Output

          if (sensorIndex .eq. 2 .and. ChannelIndex.eq.12) then 
                print*,'ln,sensorIndex, ChannelIndex',ln, SensorIndex, ChannelIndex,RTSolution(ln,m)%Brightness_Temperature
          end if

ln, sensorIndex,ChannelIndex,tb      23     2     12  257.911343594322    
ln, sensorIndex,ChannelIndex,tb      34     2     12  262.094897031404    

As you can see, for the same sensor and channel, these two methods return different TB.

In this particular test, I have total 5 instruments (simulated sections) and 1978 total channels and I selected 326 channels.

The reason I am checking this because after I passed in my channel selection into CRTM, the crtm code (ln increment after cycle channel_loop) put all selected channels in the first 326 positions in RTSolution and RTSolution_K. I thought this could potentially save memory if we only need to allocate memory for 326 channels beside improving the speed. But the computed TB is different from if my original computation (without pass in channel selection). However, If I let the ln accumulates before cycle channel_loop, the computed TB is the same as the original calculation before. It is obvious that when ln is accumulated after cycle_channel_loop, some kind of index mismatch has happened. I notice that all the calculations within the channel_loop pass in the sensorIndex and ChannelIndex which are supposed to fetch the correct crtm coefficients etc. So this is puzzling to me. I am pretty sure I passed in the correct information to CRTM since that printout is within the CRTM_K_Matrix, it got the correct ln, sensorIndex and ChannelIndex as desired.

For now, I can still pass in channel selection into CRTM and get the same values by moving the ln before the cycle channel_loop. This is a simple way to implement since the returned Solution structure matches my input measurement array. If your group has time to figure out the mysterious problem I have seen, it may potentially avoid some error if someone just use the function without evaluating the results.
"