ORAC-CC / orac

Optimal Retrieval of Aerosol and Cloud
GNU General Public License v3.0
28 stars 19 forks source link

Add flag to force nighttime retrieval. #59

Closed andyprata closed 2 years ago

andyprata commented 2 years ago

I've added a flag to force ORAC into thinking it's nighttime regardless of whether it is or it isn't. This is useful for running the retrieval during the day with only thermal channels (so far I'm finding good results when applying it to volcanic ash clouds). I force the retrieval into using nighttime cloud indexing logic by setting SPixel%Illum(:) = INight in get_indexing.F90 which subsequently triggers the cloud_indexing_logic_night routine. One potential issue is that if SPixel%Illum is written out to any of the preproc ancillary files (I haven't checked) then the flag values would be output as night when it's actually day (or twilight).

simonrp84 commented 2 years ago

Only one question: How does the user enable/disable this? It's as simple as adding Ctrl%force_nighttime_retrieval=True to the driver file, isn't it?

andyprata commented 2 years ago

Adam: Yep I agree. I did think about how to combine the do_new_night_retrieval flag with this force_nighttime_retrieval flag but didn't want to touch it until discussing with others. I think we can definitely simplify the logic that follows from case(AppCld1L) in get_indexing.F90. The solution I chose is the least invasive (i.e. everything else will work as normal unless you switch on this new flag).

The routine does consider all views at line 116-126 of get_indexing.F90 and if the illumination conditions are different then the pixel is skipped. So I assume if it's not skipped then all views have the same illumination condition and SPixel%Illum(1) can be used in the logic that follows. However, this made me realise that in the case of forcing the nighttime retrieval and setting all illumination conditions to night, we can skip this loop with:

   ! Force nighttime retrieval
   if (Ctrl%force_nighttime_retrieval) then
      SPixel%Illum(:) = INight
   else
      ! If the views have different illumination conditions then skip this pixel
      do i_view = 1, Ctrl%Ind%NViews
         SPixel%illum(i_view) = MSI_Data%illum(SPixel%Loc%X0, SPixel%Loc%Y0, i_view)
         if ((i_view > 1) .and. (.not. Ctrl%all_channels_same_view)) then
            if (SPixel%illum(i_view - 1) /= SPixel%illum(i_view)) then
               status = SPixelIndexing
               return
            end if
         end if
      end do
   end if

If you think this is OK I'll update the pull request.

Simon: Yes that's correct.

adamcpovey commented 2 years ago

Considering we literally never use anything except the first element of SPixel%illum, I think the better answer is to:

andyprata commented 2 years ago

As suggested I've now set SPixel%Illum to an integer. Main changes to the PR are:

  1. SPixel%Illum now defiend as an integer in spixel.F90 and initialised with a value of 0 in alloc_spixel.F90.
  2. SPixel%Illum array allocation/deallocation removed.
  3. MSI_Data%illum(SPixel%Loc%X0, SPixel%Loc%Y0, i_view-1) is used to compare current illumination condition to previous illumination condition in the loop of my previous comment (as SPixel%Illum is no longer an array that stores all illumination conditions).
  4. Replaced all instances of SPixel%Illum(1) with SPixel%Illum.

I've tested the changes successfully on Himawari/AHI data with force_nighttime_retrieval=True and force_nighttime_retrieval=False for the following configurations:

  1. WAT with channels 3, 4, 5, 7, 14, 15
  2. ICE with channels 3, 4, 5, 7, 14, 15
  3. EYJ with channels 3, 4, 5, 7, 14, 15
  4. EYJ with channels 13, 14, 15, 16

I've also tested the changes on some of the new LUTs and the retrieval so far works as expected.