Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

CB-TCU update #18

Merged
merged 2 commits into from
Sep 9, 2024
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
76 changes: 75 additions & 1 deletion himan-scripts/CB-TCU-cloud.lua
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,27 @@ function round(n)
return n % 1 >= 0.5 and math.ceil(n) or math.floor(n)
end

local currentProducer = configuration:GetTargetProducer()
local currentProducerName = currentProducer.GetName(currentProducer)

local filter = matrixf(9, 9, 1, missing)
--filter:Fill(1)
local kernel = {0,0,0,0,1,0,0,0,0,
0,0,1,1,1,1,1,0,0,
0,1,1,1,1,1,1,1,0,
0,1,1,1,1,1,1,1,0,
1,1,1,1,1,1,1,1,1,
0,1,1,1,1,1,1,1,0,
0,1,1,1,1,1,1,1,0,
0,0,1,1,1,1,1,0,0,
0,0,0,0,1,0,0,0,0}
filter:SetValues(kernel)

local avgkernel = {}
for i = 1, #kernel do
avgkernel[i] = kernel[i]/49
end

--Main program
--

Expand Down Expand Up @@ -59,11 +80,64 @@ if not NL or not NM or not RR then
return
end

if currentProducerName == "MEPS" or currentProducerName == "MEPSMTA" then
local Nmat = matrixf(result:GetGrid():GetNi(), result:GetGrid():GetNj(), 1, 0)
Nmat:SetValues(EL500)
EL500 = Max2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(pEL500)
pEL500 = Min2D(Nmat,filter,configuration:GetUseCuda()):GetValues()
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why take smallest pressure (highest level measured from ground up) instead of say average?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The macro is using the Max2D EL500 height above ground in the filter area. So then for the EL500 pressure I use Min2D as pressure decreases with altitude. Instead of taking the Max2D of EL500 and converting the height to pressure.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If you mean why the macro is using the Max instead of Avg I assume it is done so that airplanes stay above convection at any time. You know the tale of the horse that drowned in a river that was on average 1m deep.


Nmat:SetValues(ELmu)
ELmu = Max2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(pELmu)
pELmu = Min2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(RR)
RR = Max2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(NL)
NL = Max2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(NM)
NM = Max2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(Ttop)
Ttop = Min2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(TtopMU)
TtopMU = Min2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

filter:SetValues(avgkernel)

Nmat:SetValues(LCL500)
LCL500 = Filter2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(CAPE500)
CAPE500 = Filter2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(CIN500)
CIN500 = Filter2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(LFCmu)
LFCmu = Filter2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(pLFCmu)
pLFCmu = Filter2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(CAPEmu)
CAPEmu = Filter2D(Nmat,filter,configuration:GetUseCuda()):GetValues()

Nmat:SetValues(CINmu)
CINmu = Filter2D(Nmat,filter,configuration:GetUseCuda()):GetValues()
end

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 = 263.15 --required top T [K] (-10 degC) 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 RRlimit = 0.1 --precipitation limit [mm/h] to consider a Cb
local CAPElimit = 2.71828 --euler constant

local i = 0
Expand Down