terraref / reference-data

Coordination of Data Products and Standards for TERRA reference data
https://terraref.org
BSD 3-Clause "New" or "Revised" License
9 stars 2 forks source link

How to calculate downwelling spectral radiances #30

Closed dlebauer closed 7 years ago

dlebauer commented 8 years ago

Currently the spectral radiances are uncalibrated, and provided in the environmental logger as:

    "spectrometer": {
      "maxFixedIntensity": "16383",
      "integration time in ?s": "5000",
      "band": {
        "wavelength": 337.70483,
        "spectrum": 1500.0
      },
      "band": {
        "wavelength": 338.16013791719934,
        "spectrum": 1500.0
      },
      "band": {
        "wavelength": 338.61548740418232,
        "spectrum": 1503.0
      },
      "band": {
        "wavelength": 339.07087845402685,
        "spectrum": 1500.0
      }, ...

from 2016-04-13_00-38-15_environmentlogger.zip

Calibration files are in EnvironmentLogger/CalibrationData/ Calibrations.zip

The output of the spectrometer is 'raw' counts.

@TinoDornbusch in #26 you wrote

You need to use the attached calibration files to convert it to units of µW m-2 s-1. Careful you need to take the bandwidth of the chip into account (0.4nm) if you want to convert to µmol m-2 s-1.

Could you or @markus-radermacher-lemnatec please clarify, and provide an equation for converting the information in the calibration files to reflectances?

TinoDornbusch commented 8 years ago

The spectrometer output is irradiance. Which units do you want? W m-2 s-1?

czender commented 8 years ago

My understanding is we will multiply the calibration times the count to obtain the irradiance in uW m-2 s-1. We are interested in that, and in spectral irradiance, i.e., uW m-2 m-1 s-1 where the additional "m-1" is per unit wavelength. It's unclear to me how to obtain this. Do we divide by the bandwidth = 0.4e-9 m? That would imply that the reported irradiance only accounts for the energy in the 0.4 nm bandpass, which is fine. And with the spectral irradiance we can integrate to obtain broadband and PAR. Please confirm if this is the case and if not, precisely what to do... [NB: noticed 20160509 that I inadvertently added an extra factor of s-1 to the units in this post, and it persists in the next few messages. It should be ignored as it is folded into the W, i.e., W = J s-1. Sorry for any confusion this causes.]

TinoDornbusch commented 8 years ago

Irradiance in µW m-2 s-1 nm-1 = count * Joule/count / area / integration time

count //spectrometer output Joule/count //table in calibration file area = pi/4 *(3900µm/10^6 µm/m)^2; //this is the area of the fibre optics the light integration time = 5000µs/10^6 µs/s //spectrometer output

TinoDornbusch commented 8 years ago

Yes right......... µW m-2 s-1 nm-1 you need to divide by 0.4nm (average bandwidth).

czender commented 8 years ago

Do you have a table of the actual bandwidth per channel? Or should we effective bandwidth of each channel is 0.4 nm? We will treat these bandpasses as effective square-wave (boxcar) bandwidths unless you say otherwise.

dlebauer commented 8 years ago

Current plan is for data products to follow CF conventions, as described in #14. CF provides these variables. The units they have include steradian. Will this be required to interpret the downwelling flux?

dlebauer commented 8 years ago

@czender if I understand correctly the boxcar being integrating over [&lambda - 0.2, &lambda + 0.2]. That seems a reasonable assumption to start with - would such an assumption about the shape of this function introduce a meaningful level of error? (i.e. I don't fully understand, but the maximum flux of a boxcar would lower than a normal distribution).

@serbinsh and @ashiklom have suggested that we need FWHM for each band - is it possible that this is what 0.4um means? Does the manufacturer provide further details on the interpretation?

czender commented 8 years ago

radiance, not irradiance, has sterradian in the units. irradiance is radiance integrated over angles and multiplied by a cosine factor. my understanding is that the sensor essentially measures irradiance (at least that is what it is calibrated to return). so the docs may need to change.

i asked, and @TinoDornbusch will tell us what the 0.4 nm really means. it sounds like it is an effective squarewave bandpass, not FWHM.

TinoDornbusch commented 8 years ago

@czender Thanks for explaining the difference between radiance and irradiance.

I ask Ocean Optics for some clarification of that 0.4nm bandwidth. Nothing written in their documentation.

Note that the spectrometer has a cosine-corrected integrator. Hence we measure irradiance.

dlebauer commented 8 years ago

@czender raises the issue that the hyperspectral sensor measuring downwelling spectral radiances does not cover the range of the SWIR database.

He suggested that the key source of variability would be column water vapor (O 10%) , we could use near IR absorbance model the solar spectral flux in the longer wavebands (> 800).

My suggestion was that we could use the white reflector to calibrate an empirical model that predicts the spectral flux for all bands based on what we can observe (but is this double dipping if we also use the white reflector for calibration?). @ashiklom and @serbinsh feedback appreciated.

TinoDornbusch commented 8 years ago

@czender I asked Ocean Optics support. The bandwidth of the dispersion is NOT constant 0.4nm. To compute the width you need subtract adjacent center wavelengths.

To compute W m-2 you need to divide by this value, which roughly 0.4 nm.

dlebauer commented 8 years ago

@TinoDornbusch this is not clear to me. Could you please write out a formula? For example, if we have wavelengths

id wavelength bandwidth
1 337.70483
2 338.16013791719934
3 338.61548740418232
4 339.07087845402685

What would be the bandwidth for each of these wavelengths?

From what you say, it would be bandwidth_2 = wavelength_3 - wavelength_1 ? Or would it be (wavelength_3 - wavelength_1)/2 ? (which is closer to 0.4)?

Do you have technical documentation from Ocean Optics?

czender commented 8 years ago

@FlyingWithJerome please change the bandwith wvl_dlt from scalar constant 0.4 nm to variable bandwith as indicated above, i.e., as difference between midpoints of adjacent center wavelengths. Maybe this will solve the issue we noticed yesterday, which is that 0.4 nm gives unrealistic spectral flux peaking toward UV rather than visible wavelengths ([http://dust.ess.uci.edu/tmp/flx_spc_dwn.jpg]). I will describe this during telecon tomorrow.

TinoDornbusch commented 8 years ago

@dlebauer We have no technical documentation of Ocean Optics. I ask for it.

TinoDornbusch commented 8 years ago

@dlebauer @czender Simply making the difference between adjacent center wavelength. This came from Ocean Optics:

1 337.7048 0.45530792 2 338.1601 0.45534949 3 338.6155 0.45539105 4 339.0709 . . last pixel copy of second to last

Delta 1 = WL2-WL1 Delta 2 = WL3-WL2 . . . Last delta = second to last delta

Rounding to 3 digits is more than sufficient to achieve good accuracy.

czender commented 8 years ago

@TinoDornbusch We are having some difficulty obtaining reasonable numbers. Above you say

Irradiance in µW m-2 s-1 nm-1 = count * Joule/count / area / integration time

Are the numbers in the calibration file "joules/count", or do they include the area and integration time, so that they are "joules/count/area/integration time"? We have been assuming the latter, and so have not applied area and integration time corrections ourselves.

TinoDornbusch commented 8 years ago

The computation is correct. We have got a problem to solve. The dark reference. We need to collect a dark reference with closed sensor cap at regular intervals.

THe formula is Irradiance in µW m-2 = (count-darkref) * Joule/count / area / integration time / bandwidth

That's why your values do not make sense. At the moment we have no other solution that someone going up that gantry and covering the sensor cap at regular intervals.

I am going to ask Ocean optics how often one needs to do that dark reference. I am afraid that the chip is quite temperature sensitive and we need to do that very often. A manual solution is not applicable.

czender commented 8 years ago

@TinoDornbusch I think you neglected a factor of m-1: µW m-2 m-1 = (count-darkref) * µJoule/count / area / integration time / bandwidth Without measured spectral surface insolation, how will we calibrate the hyperspectral imagery into reflectance? The only viable paths are to "fix" the environmental logger calibration (which still won't solve the issue for SWIR calibration) or use a model instead.

czender commented 8 years ago

The broadband flux computed as I understand your instructions, is O(10^12) too big. See flx_dwn variable here. While the darkref count will reduce this, it will not eliminate it. There must be a units error or missing factor somewhere.

TinoDornbusch commented 8 years ago

I made a test myself and got 15 kW m-2 at 0:40AM. Off scale by serveal orders of magnitude. The oceanoptics Support confirmed the proposed formula. I suggest to ask Ocean optics directly for clarification:

Oliver Lischtschenko oliver.lischtschenko@oceanoptics.com

czender commented 8 years ago

I will be on travel this week and unable to pursue this further until week of 5/31. In the meantime, I encourage others to prioritize flux calibration with Ocean Optics and/or gantry system integrators.

dlebauer commented 8 years ago

@TinoDornbusch would the Ocean Optics sea breeze software help in generating calibrated output from this sensor?

dlebauer commented 8 years ago

@czender I am still waiting on a reply from Oliver but found this page on how to calibrate the sensor http://oceanoptics.com/faq/order-operations-spectral-calculations/

TinoDornbusch commented 8 years ago

I put a black cap on the spectrometer today (7.6.201) around 4PM to measure dark reference over a period of min 24h to see if there are diurnal variations. Lets hope this is quite stable since the spectrometer is in a climatized box. Let me know whether you computed values make more sense now after subtraction of dark current.

ghost commented 7 years ago

@czender - can this issue be closed?

czender commented 7 years ago

Not yet. Last time @FlyingWithJerome and I checked the irradiance was 10^12 too big. We never looked at whether subtracting the dark current measured on 20160607 improved the irradiance retrieved from the EnvironmentalLogger. We will try to look again next week. @TinoDornbusch do you now get reasonable values for irradiance when correcting for dark current?

TinoDornbusch commented 7 years ago

@czender I did not have time to make the computation. In other words, I will look at it asap.

solmazhajmohammadi commented 7 years ago

To conduct an absolute irradiance measurement, it is necessary to have the following:

Absolute irradiance (I) is computed as follows. Below, the subscript P will indicate a particular pixel for I, dL, S, D, and C. Thus, Sp refers to pixel P of the sample spectrum. Ip = (Sp - Dp) * Cp / (T * A * dLp)

Ip = (Sp - Dp) * Cp / (T * A * dLp) dL = [L(p + 1) - L(p - 1)] / 2

Rewrite W as J/sec: W/m2/nm = J/m2/sec/nm

irradiance

I checked the data with University of Oregon Solar radiation monitoring laboratory data, 0.5-2 W/m^2/nm seems to be a normal range for 400-800 nm (of course depends on when during the day we are looking at the data). Dark measurement data is missing in globus, I upload them here: https://drive.google.com/drive/folders/0Bz72iLVLzGKQQzdWcS1HY2xWT28?usp=sharing

(if the link is not opening, just copy and paste it in the url)

Please let me know if you are not getting similar results, I can check the code

czender commented 7 years ago

@FlyingWithJerome and @hmb1 pointing you to this thread. Jerome, please try to reconstruct these results with the environmental logger calculator.

FlyingWithJerome commented 7 years ago

@solmazhajmohammadi Hi Solmaz, thank you for your data and algorithm. I just downloaded these JSONs from your Google drive. Is there a variable called "dark spectrum?" Or it has been marked as something else? They seem to be the same as the older JSONs for me.

TinoDornbusch commented 7 years ago

Introduction to Environmental Biophysics says Iradiance is 500 W m-2 or 2100 µmol m-2 s-1 integrated from 400-700 nm

solmazhajmohammadi commented 7 years ago

@TinoDornbusch the graph shows the Spectral irradiance which is a function of photon wavelength. And it's true what you mentioned. Total Irradiance can be computed by integrating the spectral irradiance over all wavelengths. For the sample that I checked was 450 W/m^2 @FlyingWithJerome Dark measurment has been done on June 6th from 4pm to June 7th 10am. You can find the data here: https://drive.google.com/drive/folders/0Bz72iLVLzGKQQzdWcS1HY2xWT28?usp=sharing

solmazhajmohammadi commented 7 years ago

@FlyingWithJerome no specific name, You just dont see any changes in the spectrum (dark current is something between 1493-1507 for all the bands)

dlebauer commented 7 years ago

@solmazhajmohammadi @FlyingWithJerome is the dark spectrum something fairly constant that should be stored along with other more static calibration data (like RSR curves)? If so, could you send the file to @craig-willis so he can add it?

czender commented 7 years ago

@FlyingWithJerome since it requires major effort to perform a dark count calibration these data are likely to remain static for some time so please place directly in the environmental logger calibration routines.

solmazhajmohammadi commented 7 years ago

@dlebauer as long as we dont make any changes in the integration time, we can use the same dark measurement. I'll send the mean dark current measurement to @craig-willis today.

solmazhajmohammadi commented 7 years ago

Mean dark current values:

Spectrometer_DarkMeasurments.xlsx

solmazhajmohammadi commented 7 years ago

@czender @FlyingWithJerome

It was nice talking to you today and thank you for your feedbacks. I have plotted the raw spectrum output (unit : count) in the following figure:

count_wavelength I chose some random days and plotted them, for all of them we have a pick in the visible range.

This is a noisy sample:

spectralirradiance

You can see a spike in the plot but this is just because the collected data and dark current is noisy for the wavelength below 400nm. Here is the data from June 6 at 9:09 am: june6

Common pitfall: remember to divide out the result to the wavelength spread

dlebauer commented 7 years ago

@solmazhajmohammadi sorry I don't know all of the context, but are you saying that the time field under gantry_system_variable_metadata is provided in central time? I now realize that this is not provided according to ISO 8607 as discussed in terraref/reference-data#2 and terraref/reference-data#25; it is too late to change this field but perhaps adding a standard_time or at least a timezone field would be useful.

solmazhajmohammadi commented 7 years ago

@dlebauer Sorry, I double-checked... ,y bad It is stored at local time. I corrected the comment

czender commented 7 years ago

Thanks @solmazhajmohammadi Your plots look reasonable and have the expected spectral features. @FlyingWithJerome and I will make sure we get similar/identical results in the workflow.

FlyingWithJerome commented 7 years ago

@solmazhajmohammadi We double check the algorithms and data and yes, what we are facing now is just the plotting issue; the data themselves are OK. BTW, Did you have the result on the integration from 400 - 700nm period?

solmazhajmohammadi commented 7 years ago

For the data on June 6th, I got 310.4671 W/m^2 I plotted some sample here with the integration result using the data from 2016-07-29 spectral

ghost commented 7 years ago

@czender can this be closed?

czender commented 7 years ago

Yes, indeedy.