Developing formula for Hours of Sunshine

Very Cool! Where is this, looks like a Weather Museum!

Took it from Wikipedia. This one is in the tropics.

And here a longread how the Dutch Weather Organisation measures sunshine

http://bibliotheek.knmi.nl/knmipubWR/WR2006-06.pdf

1 Like

Dang, thatā€™s a cool interment and piece of art.

1 Like

that is a good find and an interesting read. Iā€™ll read the rest tomorrow.

A few months ago, I was looking into something very similarā€¦ I figured the instruments in the WF station could give me enough information to attempt a measurement of the clear vs. cloudiness conditions. I did a bunch of reading and tried various algorithms that used the measurements available. Originally, I looked at ratios between the UV index and solar radiation, which in theory, might give us an idea of cloud thickness. Using the premise that UV radiation is more easily scattered than the longer wavelengths in the more visible solar radiation sensorā€™s spectrum, I made a ratio between the two. The idea is that thicker clouds would scatter more UV rays so the ratio would change accordingly. Long story short, this worked some of the time but not all of the time.

I also used the expected solar radiation value and solar elevation angle values which are available to me using the ā€œastroā€ binding in OpenHAB. Google it, this is worth checking out even if you donā€™t ultimately end up using it. Itā€™s really cool!!

What I ended up using, which actually works pretty well, is to use the following variables and calculated values to come up with a the ā€œclear sky indexā€, which I found in my reading:

  • sun elevation (in degrees, from the OpenHAB astro binding) - used to calculate clear sky index
  • wind average (measured with the WF station) - used to help with identifying the fog condition
  • wind gust (measured with the WF station) - used to help with identifying fog as well
  • measured solar radiation (measured with the WF station)
  • clear day expected solar radiation (more on this below)

Essentially, what we want is to have an estimate of what the solar radiation would be if it were completely clear described as a function of solar elevation above the horizon. This is actually a pretty tricky function to get, but hereā€™s how I went about getting a version that works for my hyper-local position AND particular solar radiation sensor on the WF station.

So, as I mentioned above, I get the solar elevation angle from OpenHAB. There are other ways of getting this, potentially using the info in the links from the folks above. In order to get the expected full potential radiation at any given solar angle, there are also links above which probably can provide an estimate of what clear sky radiation would be. I tried that (also using the number derived from the OpenHAB astro binding), but found it didnā€™t match up that well to the sensor observations on a clear day. Soā€¦ what to do???

CALIBRATE THE EXPECTED RADIATION VALUE MYSELF FROM SENSOR VALUES vs. SOLAR ELEVATION ANGLE!

This only really works near the summer solstice (with data from zero to the maximum possible elevation angle at this location) on a totally clear day. I realize this might not be possible for everyone, but I found a day in early June that had just the right conditions.

As a side noteā€¦ itā€™s REALLY IMPORTANT to have the SKY device aimed as vertically as possible! I found that the temporary PVC mast I was using actually heated up in the morning and did an ā€œanti-sunflowerā€ type of movement, bending away from the sun because of the thermal expansion of the PVC on the side facing the sun. So use something that is oriented as vertically as possible and something that doesnā€™t warp when heated.

Using values from a database I am collecting, I plotted solar radiation vs. solar elevation and did a least squares fit of the data in Excel. I found that a fourth-order polynomial gave a very good fit. Hereā€™s the text of the OpenHAB ā€œruleā€ that calculates this for me. It is pretty human-readable, even if you donā€™t know the syntaxā€¦

    rule "radiation calibration"
    when
    	Item SunElevation changed
    then
    	sun_elevation 	= SunElevation.state
    	
    	if (sun_elevation > 0)
    	{
    		radiation_cal = ((c1 * Math.pow(sun_elevation.doubleValue,4)) - (c2 * Math.pow(sun_elevation.doubleValue,3)) + (c3 * Math.pow(sun_elevation.doubleValue,2)) + (c4 * sun_elevation.doubleValue))
    	}
    	else
    	{
    		radiation_cal = 0
    	}
    	postUpdate(SolarRadiationCAL, radiation_cal)
    	
    end

Ok, so now we have the ā€œradiation_calā€ value (calibrated specifically for our station and location) as a function of sun elevation angle. You then simply compare this value to the measured value from the solar sensor and generate a ratio:

measured value / radiation_cal = clear sky index ā€” or ā€œcsiā€ in the text below

Taking this ā€œcsiā€ value, along with wind values, I came up with an algorithm that tries to determine the sky conditions. Iā€™m still refining it but it actually works pretty well! Again, this is from an OpenHAB rule but the syntax should be pretty self-explanatory.

    // Cloud cover measurement algorithm
    // Todd Lorey 2018 Copyright??  :)
    // Please give me credit if you use this or publish it elewhere, but feel free to use it for yourself
    // No warranty expressed or implied for its accuracy.
    // Feel free to tweak and tune the parameters to your locality

    var Number 	csi 					= null  // clear sky index current
    var Number	csi5m 					= null  // clear sky index 5 minute average
    var Number 	csi30m 					= null  // clear sky index 30 minute average 

    var Integer	trend					= -99
    val Number	trend_threshold 		= 0.1
    var Number 	csi_delta				= null
    var Number 	csi_delta_short			= null
    var Number	wind_avg				= null
    var Number	wind_gust				= null

    val Number	clear_thresh			= 0.9
    val	Number	mostly_clear_thresh		= 0.8
    val Number	partly_cloudy_thresh 	= 0.7
    val Number	mostly_cloudy_thresh 	= 0.4
    val Number	broken_thresh			= 0.25
    val Number	overcast_thresh			= 0.05
    var	Integer	cloud_type				= -99

    var String	cloud_type_string		= "uninitialized"
    var String	trend_string			= "uninitialized"

    rule "cloud description"
    when
    	Item Weather_ClearSkyIndexLongAVG changed
    then

    csi 		= Weather_ClearSkyIndex.state as DecimalType
    csi5m 		= Weather_ClearSkyIndexAVG.state as DecimalType
    csi30m 		= Weather_ClearSkyIndexLongAVG.state as DecimalType
    wind_avg 	= Weather_WindAvg.state as DecimalType
    wind_gust	= Weather_WindGust.state as DecimalType

    csi_delta = (csi5m - csi30m)
    		if(csi_delta < 0)
    		{
    		csi_delta = ((-1) * csi_delta)
    		}
    csi_delta_short = Weather_ClearSkyIndex.maximumSince(now.minusMinutes(5), "influxdb").state as DecimalType
    	Thread::sleep(4)
    csi_delta_short = csi_delta_short - ((Weather_ClearSkyIndex.minimumSince(now.minusMinutes(5), "influxdb").state) as DecimalType)
    	Thread::sleep(4)
    	

    // Trend  //
    	trend = -98
    if ((csi5m - csi30m) > trend_threshold)
    {
    	// less cloudy
    	trend = -1
    }
    else
    {
    	if((csi30m - csi5m) > trend_threshold)
    	{
    	// more cloudy
    	trend = 1
    	}
    	else
    	{
    	// steady
    	trend = 0
    	}
    }

    // Type //
    	// Initialize //
    	cloud_type = -98
    	
    	// Fog //
    	if((0.8 > csi5m) && (csi5m > 0.15) && (0.8 > csi30m) && (csi30m > 0.15) && (trend == 0) && (wind_avg < 5) && (wind_gust < 10) && (csi_delta < 0.15) && (csi_delta_short < 0.3))
    	{
    	//Foggy
    		cloud_type = 10
    		
    		if(csi30m < 0.8) 
    		{
    		cloud_type = 11
    		}
    		if(csi30m < 0.45) 
    		{
    		cloud_type = 12
    		}
    		if(csi30m < 0.30) 
    		{
    		cloud_type = 13
    		}
    	}
    	else
    	{
    		// Clear 
    		if((csi30m >= clear_thresh) && (csi_delta_short < 0.1))
    			{
    			cloud_type = 0
    			}
    		// Mostly Clear
    		//if(((clear_thresh > csi30m) && (csi30m >= mostly_clear_thresh) && (csi_delta_short < 0.1)) || ((csi30m >= mostly_clear_thresh) && (csi_delta_short >= 0.1)))
    		if((clear_thresh > csi30m) && (csi30m >= mostly_clear_thresh) && (csi_delta_short < 0.1))
    			{
    			cloud_type = 1
    			}
    		// Partly Cloudy
    		if(((mostly_clear_thresh > csi30m) && (csi30m >= partly_cloudy_thresh) && (csi_delta_short < 0.1)) || ((csi30m >= partly_cloudy_thresh) && (csi_delta_short >= 0.1)))
    			{
    			cloud_type = 2
    			}
    		// Mostly Cloudy
    		if((partly_cloudy_thresh > csi30m) && (csi30m >= mostly_cloudy_thresh))
    			{
    			cloud_type = 3
    			}
    		// Broken
    		if((mostly_cloudy_thresh > csi30m) && (csi30m >= broken_thresh))
    			{
    			cloud_type = 4
    			}
    		// Overcast
    		if((broken_thresh > csi30m) && (csi30m >= overcast_thresh))
    			{
    			cloud_type = 5
    			}
    		// Obscured
    		if((csi30m < overcast_thresh))
    			{
    			cloud_type = 6
    			}
    		
    	}

    switch(cloud_type) 
    		
    	{
    	case 10 :
    		{
    		cloud_type_string = "Foggy"
    		}
    	case 11 :
    		{
    		cloud_type_string = "Thin fog"
    		}
    	case 12 :
    		{
    		cloud_type_string = "Medium fog"
    		}
    	case 13 :
    		{
    		cloud_type_string = "Heavy fog"
    		}
    	case 0 :
    		{
    		cloud_type_string = "Clear"
    		}
    	case 1 :
    		{
    		cloud_type_string = "Mostly clear"
    		}
    	case 2 :
    		{
    		cloud_type_string = "Partly cloudy"
    		}
    	case 3 :
    		{
    		cloud_type_string = "Mostly cloudy"
    		}
    	case 4 :
    		{
    		cloud_type_string = "Broken clouds"
    		}
    	case 5 :
    		{
    		cloud_type_string = "Overcast"
    		}
    	case 6 :
    		{
    		cloud_type_string = "Obscured"
    		}
    	case -98 :
    		{
    		cloud_type_string = "UNKNOWN"
    		}
    	case -99 :
    		{
    		cloud_type_string = "UNINITIALIZED"
    		}
    	}

    switch(trend) 
    		
    	{
    	case -1 :
    		{
    		trend_string = "decreasing coverage"
    		}
    	case 0 :
    		{
    		trend_string = "steady coverage"
    		}
    	case 1 :
    		{
    		trend_string = "increasing coverage"
    		}
    	case -98 :
    		{
    		trend_string = "UNKNOWN"
    		}
    	case -99 :
    		{
    		trend_string = "UNINITIALIZED"
    		}
    	}


    if((cloud_type == 0) && (trend == 0))
    {
    	postUpdate(Weather_CloudString, cloud_type_string)
    	//logInfo("info", csi.toString+ ", " +csi5m.toString+ ", " +csi30m.toString+ ", " +csi_delta.toString+ ", " +csi_delta_short.toString+ ", " +cloud_type_string)
    }
    else
    {
    	postUpdate(Weather_CloudString, cloud_type_string+ ", " +trend_string)
    	//logInfo("info", csi.toString+ ", " +csi5m.toString+ ", " +csi30m.toString+ ", " +csi_delta.toString+ ", " +csi_delta_short.toString+ ", " +cloud_type_string+ ", " +trend_string )
    }
    end	

I hope this helps someone out if they are trying to do something similar. I canā€™t promise much explanation or support in implementing this for yourself. It does require some knowledge of programming and lots of head scratching to get it to work.

6 Likes

Excellent find. Still reading, but if:

/*
Further investigations resulted in a mean value for the threshold of 120 Wm-2, which was accepted in 1989 by the WMO as the actual threshold. As a reference sensor for the detection of the threshold irradiance, a pyrheliometer was recommended, an instrument that measures the direct normal solar irradiance (DNSI).
*/

The above document then goes on to talk about the Slob/Bergman algorythm, and how itā€™s used to massage pyranometer readings into the above.

Neat!

Currently, @tlwolter and I have a base UV index range that works pretty well at ā€˜predictingā€™ the current sun/cloud cover.

Looking at implementing it as a gaussian distribution, combined with NOAA UV monthly indexes for the US, it looks like we should be able to create an algorithm to predict the current conditions pretty well. Combining that with current solar radiation numbers should give us a good idea on the ā€˜hours of sunshineā€™ each day.

Gonna repost with test results soon.

ā€“Sam

2 Likes

i donā€™t understand why this is moved it to developers, as it was really meant to be a feature request. The total hours of sunshine each day, is nice info for casual weather geeks. It would be great if it was just one of the derived values that is shown in the app.

of course it is fine if some third party software would also show it, but I really would like it to have it in the existing app. @dsj

I asked that it be split but all the posts were moved.

We should not have hijacked your Feature Request post.

@sunny,

To clarify, are you wanting the total two the the sky is clear and the sun shines or the total daylight hours?

ā€“Sam

Total daylight hours is pretty simple.

I think we should work on Total hours of sunshine.
Total hours of cloud overhead sounds interesting but could be difficult.

Sorry if my geeking out caused this thread to be derailed. I meant no harm, only wanted to share some thoughts on a related subject.

3 Likes

Itā€™s not an issue Todd. Thatā€™s why I asked that it be split. The OP can continue to advocate for this feature from WeatherFlow and we can work on a formula.

Hi,

Personally, I just like to have the total hours the sun is shining (derived from the measured solar radiance). I did not ask for the hours of daylight (which is an astronomical value based only on location and time of the year).

So it the hours one could actually see the sun and it is not obscured by the clouds.

In the beginning of the thread someone mentioned a collection of javascript that calculates the hours of daylight. I assume it is a language thing as in my language ā€˜hours of sunshineā€™ would be very clear.

However, as mentioned, it might be true that in order to derive the hours of sunshine from the measured.radiance, other values are needed like, time of the day, time of the year, location, For example, close to sunset, the sun might be shining, but the solar radiance values might be less then on a bright day at noon but the sun being blocked by a small cloud and it isnā€™t shining. Or, as someone suggested, it might be derived from solar radiance and uv radiance.

2 Likes

I originally misunderstood what you wanted. The hours between sunrise and sunset are easy. Niw we need to figure out the number of minutes that the sun actually reaches the sensor on the Sky, or ā€œHow much time did the sun shine on the ground?ā€

We just need to agree on minimum levels of UV, illumination and solar radiation.

No problem, Todd, I like the discussion about the best formula to calculate the hours of sunshine. I just didnā€™t think it needed splitting or the accidental complete move of the thread. But if people think that finding the best formula is better handled in a different thread, fine. I hope the end result will be moved or copied back to the feature request, so weatherflow staff knows what to implement to resolve this request. After all, a feature request that includes the right formula is probably a lot stronger.

1 Like

One option would be : open a new thread in ā€œfeature requestsā€ and Linke it to this thread.
What do you think ?

Iā€™m not sure if we can use the difference between uv index and solar radiance. Logically I would say yes, but what I have read so far, using spectral differences isnā€™t something that is being done in the instruments that replaced the Campbell-stokes instruments. There they donā€™t use uv index.

1 Like

Letā€™s just wait and see where this thread ends and then we can copy the result back. Itā€™s already confusing enough.

1 Like