Skip to content

Commit

Permalink
Merge branch 'master' of github.com:fmidev/himan
Browse files Browse the repository at this point in the history
  • Loading branch information
@mpartio
mpartio committed May 7, 2024
2 parents 090008e + af90d90 commit 1fbbf59
Showing 1 changed file with 56 additions and 44 deletions.
100 changes: 56 additions & 44 deletions himan-scripts/CB-TCU-cloud.lua
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,31 +10,39 @@ local MU = level(HPLevelType.kMaximumThetaE,0)
local HL = level(HPLevelType.kHeightLayer,500,0)
local HG = level(HPLevelType.kHeight,0)

EL500 = luatool:Fetch(current_time, HL, param("EL-LAST-HPA"), current_forecast_type)
LCL500 = luatool:Fetch(current_time, HL, param("LCL-HPA"), current_forecast_type)
LFC500 = luatool:Fetch(current_time, HL, param("LFC-HPA"), current_forecast_type)
CAPE500 = luatool:Fetch(current_time, HL, param("CAPE-JKG"), current_forecast_type)
CIN500 = luatool:Fetch(current_time, HL, param("CIN-JKG"), current_forecast_type)

LCLmu = luatool:Fetch(current_time, MU, param("LCL-HPA"), current_forecast_type)
LFCmu = luatool:Fetch(current_time, MU, param("LFC-HPA"), current_forecast_type)
ELmu = luatool:Fetch(current_time, MU, param("EL-LAST-HPA"), current_forecast_type)
CINmu = luatool:Fetch(current_time, MU, param("CIN-JKG"), current_forecast_type)
CAPEmu = luatool:Fetch(current_time, MU, param("CAPE-JKG"), current_forecast_type)

Ttop = luatool:Fetch(current_time, HL, param("EL-K"), current_forecast_type)
Tbase = luatool:Fetch(current_time, HL, param("LCL-K"), current_forecast_type)
TtopMU = luatool:Fetch(current_time, MU, param("EL-K"), current_forecast_type)
--LFC probably better than LCL for elev conv. base
TbaseMU = luatool:Fetch(current_time, MU, param("LFC-K"), current_forecast_type)

if not EL500 or not LCLmu or not Ttop then
local EL500 = luatool:Fetch(current_time, HL, param("EL-LAST-M"), current_forecast_type)
local pEL500 = luatool:Fetch(current_time, HL, param("EL-LAST-HPA"), current_forecast_type)
local LCL500 = luatool:Fetch(current_time, HL, param("LCL-M"), current_forecast_type)
local CAPE500 = luatool:Fetch(current_time, HL, param("CAPE-JKG"), current_forecast_type)
local CIN500 = luatool:Fetch(current_time, HL, param("CIN-JKG"), current_forecast_type)
local LFCmu = luatool:Fetch(current_time, MU, param("LFC-M"), current_forecast_type)
local pLFCmu = luatool:Fetch(current_time, MU, param("LFC-HPA"), current_forecast_type)
local ELmu = luatool:Fetch(current_time, MU, param("EL-LAST-M"), current_forecast_type)
local pELmu = luatool:Fetch(current_time, MU, param("EL-LAST-HPA"), current_forecast_type)
local CINmu = luatool:Fetch(current_time, MU, param("CIN-JKG"), current_forecast_type)
local CAPEmu = luatool:Fetch(current_time, MU, param("CAPE-JKG"), current_forecast_type)
local Ttop = luatool:Fetch(current_time, HL, param("EL-K"), current_forecast_type)
local TtopMU = luatool:Fetch(current_time, MU, param("EL-K"), current_forecast_type)

if not EL500 or
not pEL500 or
not LCL500 or
not CAPE500 or
not CIN500 or
not LFCmu or
not pLFCmu or
not ELmu or
not pELmu or
not CINmu or
not CAPEmu or
not Ttop or
not TtopMU then
logger:Error("Some data not found")
return
end

NL = luatool:Fetch(current_time, HG, param("NL-PRCNT"), current_forecast_type)
NM = luatool:Fetch(current_time, HG, param("NM-PRCNT"), current_forecast_type)
local NL = luatool:Fetch(current_time, HG, param("NL-PRCNT"), current_forecast_type)
local NM = luatool:Fetch(current_time, HG, param("NM-PRCNT"), current_forecast_type)

if not NL then
NL = luatool:Fetch(current_time, HG, param("NL-0TO1"), current_forecast_type)
Expand All @@ -44,16 +52,19 @@ if not NM then
NM = luatool:Fetch(current_time, HG, param("NM-0TO1"), current_forecast_type)
end

RR = luatool:Fetch(current_time, HG, param("RRR-KGM2"), current_forecast_type)
local RR = luatool:Fetch(current_time, HG, param("RRR-KGM2"), current_forecast_type)

CBlimit = 12 --required vertical thickness [degrees C] to consider a CB (tweak this..!)
TCUlimit = 9 --required vertical thickness [degrees C] to consider a TCU (tweak this..!)
CBtopLim = 263.15 --required top T [K] to consider a CB (tweakable!)
CINlimTCU = -1 --CIN limit for TCu
RRlimit = 0.1 -- precipitation limit [mm/h] to consider a Cb
if not NL or not NM or not RR then
logger:Error("Some data not found")
return
end

--Max height [FL] to check for top
TopLim = 650
local CBlimit = 2000 --required vertical thickness [m] to consider a CB (tweak this..!)
local TCUlimit = 1500 --required vertical thickness [m] to consider a TCU (tweak this..!)
local CBtopLim = -10 --required top T [K] to consider a CB (tweakable!)
local CINlimTCU = -1 --CIN limit for TCu
local RRlimit = 0.1 -- precipitation limit [mm/h] to consider a Cb
local CAPElimit = 2.71828 --euler constant

local i = 0
local res = {}
Expand All @@ -64,10 +75,11 @@ for i=1, #EL500 do
res[i] = Missing

--TCU
if ((Tbase[i]-Ttop[i]>TCUlimit) and (NL[i]>0) and (CIN500[i]>CINlimTCU)) then
res[i] = FlightLevel_(EL500[i]*100)
if ((EL500[i] - LCL500[i] > TCUlimit) and (NL[i] > 0) and (CIN500[i] > CINlimTCU) ) then
--we don't use vertical search for flight level of EL500 but calculate directly from EL500 pressure
res[i] = FlightLevel_(pEL500[i] * 100)
--Limit top value
if (CAPE500[i]>math.exp(1)) then
if (CAPE500[i] > CAPElimit) then
--Add for overshooting top based on CAPE
res[i] = -(res[i] + CAPE500[i] / (math.log(CAPE500[i]) * 10))
else
Expand All @@ -76,33 +88,33 @@ for i=1, #EL500 do
end

--CB
if ((Ttop[i]<CBtopLim) and (Tbase[i]-Ttop[i]>CBlimit) and (RR[i]>RRlimit)) then
res[i] = FlightLevel_(EL500[i]*100)
if ((Ttop[i] < CBtopLim) and (EL500[i] - LCL500[i] > Cblimit) and (RR[i] > RRlimit)) then
res[i] = FlightLevel_(pEL500[i] * 100)
--Limit top value
if (CAPE500[i]>math.exp(1)) then
if (CAPE500[i] > math.exp(1)) then
--Add for overshooting top based on CAPE
res[i] = res[i] + CAPE500[i] / (math.log(CAPE500[i]) * 10)
end
end

--If no TOP from above, check also with MU values, for elev. conv. only from blw 3,5km
if ( IsMissing(res[i]) and (TbaseMU[i]<Tbase[i]) and (LFCmu[i]>650)) then
-- TCU
if ((TbaseMU[i]-TtopMU[i]>TCUlimit) and ((NL[i]>0) or (NM[i]>0)) and (CINmu[i]>CINlimTCU)) then
res[i] = FlightLevel_(ELmu[i]*100)
if ( IsMissing(res[i]) and (LFCmu[i] > LCL500[i]) and (pLFCmu[i] > 650)) then
-- TCU elevated
if ((ELmu[i] - LFCmu[i] > TCUlimit) and ((NL[i] > 0) or (NM[i] > 0)) and (CINmu[i] > CINlimTCU)) then
res[i] = FlightLevel_(pELmu[i] * 100)
--Limit top value
if (CAPEmu[i]>math.exp(1)) then
if (CAPEmu[i] > CAPElimit) then
--Add for overshooting top based on CAPE, +1000ft/350J/kg (tweak this!)
res[i] = -(res[i] + CAPEmu[i] / (math.log(CAPEmu[i]) * 10))
else
res[i] = -res[i]
end
end
--CB
if ((TtopMU[i]<CBtopLim) and (TbaseMU[i]-TtopMU[i]>CBlimit) and (RR[i]>RRlimit)) then
res[i] = FlightLevel_(ELmu[i]*100)
--CB elevated
if ((TtopMU[i] < CBtopLim) and (ELmu[i] - LFCmu[i] > CBlimit) and (RR[i] > RRlimit)) then
res[i] = FlightLevel_(pELmu[i] * 100)
--Limit top value
if (CAPEmu[i]>math.exp(1)) then
if (CAPEmu[i] > CAPElimit) then
--Add for overshooting top based on CAPE, +1000ft/350J/kg (tweak this!)
res[i] = res[i] + CAPEmu[i] / (math.log(CAPEmu[i]) * 10)
end
Expand Down

0 comments on commit 1fbbf59

Please sign in to comment.