nasa / HyperCP

Other
41 stars 23 forks source link

NaN and 0 values returned by 6S #205

Open raphidoc opened 4 months ago

raphidoc commented 4 months ago

At some wavelength, 6S might return NaN or 0 values depending on the output quantity. For example, the percentage of direct and diffuse irradiance at a certain wavelength would be returned as 0 if the absolute irradiance at the same wavelength is returned as NaN.

The problem is that 0 might be an appropriate value for a physical quantity. For example, environmental irradiance is often returned as 0 over a large part of the spectrum. In that case, it effectively indicates a 0 contribution from this irradiance component. If only the various ratios are reported as 0 instead of NaN, this might help to distinguish them.

The other point is that 6S run on 20 node wavelength and then perform interpolation of the values at the requested wavelength. The interpolation formula is reported as: image

The 20 node wavelength are:

      data wldisc /0.350,0.400,0.412,0.443,0.470,0.488,0.515,0.550,
     s             0.590,0.633,0.670,0.694,0.760,0.860,1.240,1.536,
     s             1.650,1.950,2.250,3.750/

I wonder if it would be possible to run 6S on those 20 wavelengths only and then interpolate to the instrument spectrum?

oceancolorcoder commented 4 months ago

If irradiance is 0, then I would think direct/diffuse could go to NaN, but why would absolute irradiance return NaN at some wavelengths instead of 0? What distinguishes irradiance of 0 and NaN in 6S?

I don't quite follow your second question. Based on your question, it seems like 6S is already run on only those 20 wavelengths and then interpolated to the instrument. How is what you are suggesting different?

raphidoc commented 4 months ago

Regarding the first question, I don't know why 6S returns NaN at certain wavelengths. I will have to run some tests to better understand this issue and the Fortran code deals with NaN and 0.

I'll try to clarify my second question:

  1. 6S takes a wavelength parameter as input, according to the 6S code documentation, it can take the following forms:

        iwave input of the spectral conditions
    
        you choose to define your own spectral conditions: iwave=-1,0 or 1
        (three user s conditions )
    
        -2  enter wlinf, wlsup, the filter function will be equal to 1
        over the whole band (as iwave=0) but step by step output will be printed
    
        -1  enter wl (monochr. cond,  gaseous absorption is included)
    
        0  enter wlinf, wlsup. the filter function will be equal to 1 over the whole band.
    
        1  enter wlinf, wlsup and user's filter function s(lambda) ( by step of 0.0025 micrometer).
  2. The -1 condition is currently in use by HyperCP with the call to Py6S.SixSHelpers.Wavelengths.run_wavelengths(s, 1e-3*wvl, n=n_cores)

  3. So 6S is run one time for each wvl. For each requested wavelength, 6S runs the model on 20 wavelengths and then interpolates to the requested wavelength.

With my approach, running 6S 255 times takes 0.4 seconds, and running it 20 times takes 0.04 seconds. If we could perform the same interpolation as within 6S from the 20 node wavelengths, I expect that this would be faster than running 6S 255 times. It also could solve the problem of 6S returning NaN at certain wavelengths, but again I don't know why this happen in the first place.

oceancolorcoder commented 3 months ago

@raphidoc I see now. I agree that it would make far more sense to run 6S once on the native 20 bands and interpolate those to our hyperspectral data rather than running 6S for each of our hyperspectral bands. We will have to see which (if any) of the native 20 bands are returned with NaNs to decide what to do about it.

oceancolorcoder commented 3 months ago

NaN catch and patch for 6S results addressed in PR #210

Leaving issue open until we can discuss the 6S interpolation with the group next week.

raphidoc commented 3 months ago

I just noticed something: three wavelengths at which NaN are returned are: 443.05, 516.13, 589.8, 633.28, 633.45, 670.01, 670.13

Those share the same integer part as three of the node wavelengths of 6S

data wldisc /0.350,0.400,0.412,0.443,0.470,0.488,0.515,0.550,
     s             0.590,0.633,0.670,0.694,0.760,0.860,1.240,1.536,
     s             1.650,1.950,2.250,3.750/
oceancolorcoder commented 3 months ago

Weird. Same result with your py6S_json wrapper?

raphidoc commented 3 months ago

Yes, this comes from 6S itself.

oceancolorcoder commented 3 months ago

Perhaps a long shot, but is the 6S team responsive to feedback and issues submission? Can you test another band (e.g., spoof the input wavelength) to confirm the weirdness?

raphidoc commented 3 months ago

The three wavelengths reported above trigger the bug, but it doesn't seem related to the number or values of the digits. I tested with various wavelengths, but I couldn't see a pattern. All the 20 native wavelengths of 6S produce a correct output.

I wrote a mail to Svetlana Y. Kotchenova, one of the authors of the software.