Closed fherreazcue closed 5 months ago
The diagnostic names that can in principle be extracted from rip
This is the code from the NCL routine to get the height data:
if( any( variable .eq. (/"geopt","geopotential","z","height", "height_agl"/) ) ) then
; Height [=full geopotentail height / 9.81]
if(isfilevar(nc_file,"PH")) then
var = _get_wrf_var(file_handle,"PH",time)
PHB = _get_wrf_var(file_handle,"PHB",time)
var = var + PHB
z = wrf_user_unstagger(var,var@stagger)
z@description = "Geopotential"
elseif(isfilevar(nc_file,"GHT")) then
;; may be a met_em file - see if we can get GHT - data in met_em file is Height in M
z = _get_wrf_var(file_handle,"GHT",time)
z = z * 9.81
z@description = "Geopotential"
z@units = "m2 s-2"
else
print("wrf_user_getvar: Error: cannot find 'PH' or 'GHT' on the file")
return
end if
if( any( variable .eq. (/"z","height", "height_agl"/) ) ) then
z = z / 9.81
z@description = "Height"
if (variable .eq. "height_agl") then
if (isfilevar(nc_file,"HGT") ) then
ter = _get_wrf_var(file_handle, "HGT", time)
rank := dimsizes(dimsizes(ter))
if (rank .eq. 2) then
ter_conform = conform_dims(dimsizes(z), ter, (/1,2/))
elseif (rank .eq. 3) then
ter_conform = conform_dims(dimsizes(z), ter, (/0,2,3/))
elseif (rank .eq. 4) then
ter_conform = conform_dims(dimsizes(z), ter, (/0,1,3,4/))
elseif (rank .eq. 5) then
ter_conform = conform_dims(dimsizes(z), ter, (/0,1,2,4,5/))
elseif (rank .eq. 6) then
ter_conform = conform_dims(dimsizes(z), ter, (/0,1,2,3,5,6/))
end if
z = z - ter_conform
z@description = "Height_AGL"
else
print("wrf_user_getvar: Error: cannot find 'HGT' in the file")
return
end if
end if
z@units = "m"
end if
return(z)
end if
Looks like the process to get height above sea level is: 1) Add PH + PHB 2) Unstagger the grid of the resultant data 3) Divide by gravitational acceleration rate (9.81 m/s^2)
Ahh - looks like RIP already returns GHT though - so perhaps just use that.
(I had assumed these were wrfout variables listed here, but perhaps RIP does this postprocessing work too)
From Dave: """ Can you output both ght and ghtagl? By mixing ratio, I implied water vapor mixing ratio, so variable qvp. Rather than tdk and tmk, could you do tdp and tmc? I forgot about all those variables. As long as we're doing that, could you output these, too? lcl, lfc, omg, pvm, pvo, sateth, stb, stbe, stbz, tdd, wsp, wdr, www Looks like you'll be able to get CIN and CAPE easily (cin3 and cap3, mcap and mcin). """
All diagnostics added in https://github.com/UoMResearchIT/SpanishPlumeAnalysis/commit/beef8cfefe374e55bbb6f63ac51469313e6ed38a, ready to be run on csf. However, some issues with memory were encountered on local tests (see this comment), so it might be that we need to run the trajectory calculation twice to get all the desired diagnostics (running the trajectory with only the second half of diagnostics does work, so it must be that it runs out of memory). Before changing the code though, I need to see if the same problem is encountered when ran on the csf.
The issue persisted on the csf, so it was necessary to split the run in two (see https://github.com/UoMResearchIT/SpanishPlumeAnalysis/commit/7d0572a345cab453a5da91e6b9d33c80563f039c and https://github.com/UoMResearchIT/SpanishPlumeAnalysis/commit/2d6adab34a3be94b18d3cbd5b8187607288e88d8). Trajectories have been generated, with outputs for csvs split in files with g1 or g2 in the name, and the results are in dropbox.
CSV files along trajectories of 90h ran successfully, results are in Dropbox.
Backward trajectories from Nottingham at
Ending locations at 950, 925, 900, 850, 800, 750, 700, 650, 600, 500, 400, 300, and 250 hPa.
Simulations: Control, ALB-40, ALB-Std, ALB-90
Plots and CSV file with following columns: