--Variables à éditer--------------
local Ville = "pws:ILAGACIL2" --Ville de référence
local idxLux = '320' --Capteur virtuel de Luminosité
local idxOcta = '5' -- idx de la variable Utilisateur
local wuAPIkey = "48a08328a93a18a1"
local DEBUG = 'OUI'

--import des fontions pour lire le JSON
json = (loadfile "/opt/domoticz/scripts/lua/JSON.lua")()

-- Fonction déterminant si année bissextile
function AnneeBissextile(annee)
   return annee%4==0 and (annee%100~=0 or annee%400==0)
end

---------------------------------
-- Fonction spliter
function split(s, delimiter)
   result = {};
   for match in (s..delimiter):gmatch("(.-)"..delimiter) do
      table.insert(result, match);
   end
   return result;
end


---------------------------------
-- Fonction Arrondir
function arrondir(num, dec)
   if num == 0 then
      return 0
   else
      local mult = 10^(dec or 0)
      return math.floor(num * mult + 0.5) / mult
   end
end

---------------------------------
-- Procedure principale
---------------------------------

---------------------------------
-- Meteo API Wunderground
local config=assert(io.popen('curl http://api.wunderground.com/api/'..wuAPIkey..'/conditions/q/FR/'..Ville..'.json'))
local location = config:read('*all')
config:close()
local jsonLocation = json:decode(location)

local Latitude = jsonLocation.current_observation.display_location.latitude
local Longitude = jsonLocation.current_observation.display_location.longitude
local Altitude = jsonLocation.current_observation.display_location.elevation
PressionRelative = jsonLocation.current_observation.pressure_mb

----------------------------------
local An = os.date("%Y")
local NiemeJourDeLAnnee = os.date("%j")
---------------------------------
DateHeure = os.date("%Y-%m-%d %H:%M:%S", os.time())
-- Début debug
if( DEBUG == 'OUI') then
   print("=====================================")
   print(os.date("%Y-%m-%d %H:%M:%S", os.time()))
   print(Ville .. ", " .. Latitude .. ", " .. Longitude)
   print("Altitude = " .. tostring(Altitude) .. " m")
   print("NiemeJourDeLAnnee = " .. NiemeJourDeLAnnee)
end
---------------------------------
--calcul de l'année bissextile
if  AnneeBissextile(An) == true then
   if( DEBUG == 'OUI') then
      print( An .. " est bissextile.")
   end
   JourDansLAnnee = 366
else
   if( DEBUG == 'OUI') then
      print( An .. " n'est pas bissextile.")
   end
   JourDansLAnnee = 365
end
---------------------------------
-- Vitesse angulaire = Combien de degrés par jour
VitesseAngulaire = 360/365.25 ----JourDansLAnnee -- ou approximativement 365.25
if( DEBUG == 'OUI') then
   print("Vitesse angulaire = " .. VitesseAngulaire .. " par jour")
end
---------------------------------
-- Formule Declinaison = ArcSin(0,3978 x Sin(Va x (j - (81 - 2 x Sin(Va� x (j - 2))))))
local Declinaison = math.deg(math.asin(0.3978 * math.sin(math.rad(VitesseAngulaire) *(NiemeJourDeLAnnee - (81 - 2 * math.sin((math.rad(VitesseAngulaire) * (NiemeJourDeLAnnee - 2))))))))
if( DEBUG == 'OUI') then
   print("La déclinaison = " .. Declinaison .. "°")
end
---------------------------------
-- Temps universel décimal (UTC)
TempsDecimal = (os.date("!%H") + os.date("!%M") / 60)
if( DEBUG == 'OUI') then
   print("Temps universel decimal (UTC)".. TempsDecimal .." H.dd")
end
---------------------------------
-- Temps solaire
HeureSolaire = TempsDecimal + (4 * Longitude / 60 )
if( DEBUG == 'OUI') then
   print("Temps solaire ".. HeureSolaire .." H.dd")
end
---------------------------------
-- Angle horaire du soleil
AngleHoraire = 15 * ( 12 - HeureSolaire )
if( DEBUG == 'OUI') then
   print("Angle Horaire = ".. AngleHoraire .. "°")
end
---------------------------------
-- La hauteur du soleil (Elévation ou altitude)
HauteurSoleil = math.deg(math.asin(math.sin(math.rad(Latitude))* math.sin(math.rad(Declinaison)) + math.cos(math.rad(Latitude)) * math.cos(math.rad(Declinaison)) * math.cos(math.rad(AngleHoraire))))
if( DEBUG == 'OUI') then
   print("Hauteur du soleil = " .. HauteurSoleil .. "°")
end
local Azimut = math.acos((math.sin(math.rad(Declinaison)) - math.sin(math.rad(Latitude)) * math.sin(math.rad(HauteurSoleil))) / (math.cos(math.rad(Latitude)) * math.cos(math.rad(HauteurSoleil) ))) * 180 / math.pi
local SinAzimut = (math.cos(math.rad(Declinaison)) * math.sin(math.rad(AngleHoraire))) / math.cos(math.rad(HauteurSoleil))
if(SinAzimut<0) then
   Azimut=360-Azimut
end
if( DEBUG == 'OUI') then
   print("Azimut du soleil = " .. Azimut .. "°")
end
---------------------------------
-- La durée d'insolation journalière - non stockée en VG
DureeInsolation = math.deg(2/15 * math.acos(- math.tan(math.rad(Latitude)) * math.tan(math.rad(Declinaison))))
DureeInsolation = arrondir(DureeInsolation,2)
if( DEBUG == 'OUI') then
   print("La durée d'insolation journalière = " .. DureeInsolation .." H.dd")
end
---------------------------------
-- Constantes Solaire
ConstanteRatiationSolaire = 1361 -- W/m²
ConstanteRadiationLux = 200000 -- Lux
---------------------------------
-- Rayonnement solaire (en W/m²) présent à l'entrée de l'atmosphère.
RadiationAtm = ConstanteRatiationSolaire * (1 +0.034 * math.cos( math.rad( 360 * NiemeJourDeLAnnee / JourDansLAnnee )))
if( DEBUG == 'OUI') then
   print("Radiation max en atmosphère = " .. arrondir(RadiationAtm,2) .. " W/m²")
end
---------------------------------
-- Coefficient d'attenuation M
PressionAbsolue = PressionRelative - arrondir((Altitude/ 8.3),1) -- hPa
SinusHauteurSoleil = math.sin(math.rad(HauteurSoleil))
M0 = math.sqrt(1229 + math.pow(614 * SinusHauteurSoleil,2)) - 614 * SinusHauteurSoleil
M = M0 * PressionRelative/PressionAbsolue
if( DEBUG == 'OUI') then
   print("Pression relative locale = " .. PressionRelative .. " hPa")
   print("Pression absolue atmosphère = " .. PressionAbsolue .. " hPa")
   print("Coefficient d'attenuation = " .. M )
end
---------------------------------
-- Récupérer message SYNOP sur le site Ogimet
heureUTCmoins1 = os.date("!%H")-1
if string.len(heureUTCmoins1) == 1 then
   heureUTCmoins1 = "0" .. heureUTCmoins1
end
UTC = os.date("%Y%m%d").. heureUTCmoins1.."00" -- os.date("!%M")
local WMOID = jsonLocation.current_observation.display_location.wmo
print("WMOID="..WMOID)
local ogimet=assert(io.popen('curl "http://www.ogimet.com/cgi-bin/getsynop?block='..WMOID..'&begin='..UTC..'"'))
local synop = ogimet:read('*all')
ogimet:close()
rslt = split(synop,",")
CodeStation = rslt[1]
rslt = split(synop, " "..CodeStation.. " ")
Trame = string.gsub(rslt[2], "=", "")
Trame = CodeStation .." ".. Trame
rslt = split(Trame, " ")
if( DEBUG == 'OUI') then
   print("Horodatage UTC = " .. UTC)
   print("Station SYNOP = " .. WMOID)
end
---------------------------------
-- Récupérer le premier caractere du 3eme mot = Nebulosité en Octa
Octa = string.sub(rslt[3], 1, 1)
--print(rslt[3])
-- 0         Pas de couverture nuageuse
-- 1-8       Huitième
-- 9         Brouillard
-- /        Couverture indiscernable
-- cas particulier si valeur indéterminé un slash est renvoyé.
if Octa == "/" then
   -- si la couverture est indiscernable, on reprend l'ancienne valeur
   -- cela evite les violentes cassures dans la courbe
   --Octa = uservariables['octa']
   Octa = arg[1]
elseif Octa == "9" then
   Octa = 8
end
if( DEBUG == 'OUI') then
   print( Octa .. " Octa")
end
-- stockage de la variable octa
--commandArray['Variable:octa']=tostring(Octa)
os.execute('/usr/bin/curl "http://127.0.0.1:8080/json.htm?type=command&param=updateuservariable&idx='..idxOcta..'&vname=octa&vtype=0&vvalue='..tostring(Octa)..'"')


---------------------------------
-- Facteur d'atténuation des couches nuageuses Kc
-- Kc=1-(0.75*((OCTA)**(3.4))
Kc=1-0.75*(math.pow(Octa/8,3.4))
if( DEBUG == 'OUI') then
   print("Kc = " .. Kc)
end
---------------------------------
-- Au lever/coucher du soleil, on atteind les limites de précisions de ces calculs.
-- J'interrompts donc le calcul de radiation dès 1°.
if HauteurSoleil > 1 then
-- Radiation directe
   RadiationDirecte = RadiationAtm * math.pow(0.6,M) * SinusHauteurSoleil
-- Radiation Diffuse
   RadiationDiffuse = RadiationAtm * (0.271 - 0.294 * math.pow(0.6,M)) * SinusHauteurSoleil
-- Radiation totale
   RadiationTotale = RadiationDiffuse + RadiationDirecte
-- Radiation en Lux : --  1 Lux = 0,0079 W/m²
   Lux = RadiationTotale / 0.0079
   --Lux = ConstanteRadiationLux / ConstanteRatiationSolaire * RadiationTotale
-- Le rayonnement solaire avec ciel nuageux
   RTOTC = RadiationTotale * Kc
-- Radiation en Lux pondéré
   -- LuxPondere = ConstanteRadiationLux / ConstanteRatiationSolaire * RTOTC
   LuxPondere = RTOTC / 0.0079
   print("Radiation totale en lux pondéré = ".. arrondir(LuxPondere,2).." Lux")
   if( DEBUG == 'OUI') then
      print("RadiationDirecte = ".. arrondir(RadiationDirecte,2) .." W/m²")
      print("Radiation Diffuse = ".. arrondir(RadiationDiffuse,2) .." W/m²")
      print("Radiation totale = " .. arrondir(RadiationTotale,2) .." W/m²")
      print("Radiation totale en lux = ".. arrondir(Lux,2).." Lux")
      print("Le rayonnement solaire avec pondération = " .. arrondir(RTOTC,2))
   end
else
   --le soleil est trop bas
   RTOTC = 0
   LuxPondere = 0
end
   --mise à jour du dispositif
   os.execute('/usr/bin/curl "http://127.0.0.1:8080/json.htm?type=command&param=udevice&idx='..idxLux..'&svalue="'..tostring(LuxPondere)..' -s -o /dev/null')