MITgcmTools.jl
Set of tools for the analysis, manipulation, etc of MITgcm input and output files.
MITgcmTools.MatrixInterp
MITgcmTools.MetaFileRead
MITgcmTools.MixedLayerDepth
MITgcmTools.SeaWaterDensity
MITgcmTools.convert2array
MITgcmTools.convert2gcmfaces
MITgcmTools.findtiles
MITgcmTools.parsemeta
MITgcmTools.prep_MTRX
MITgcmTools.qwckplot
MITgcmTools.qwckplot
MITgcmTools.qwckplot
MITgcmTools.readAvailDiagnosticsLog
MITgcmTools.read_SPM
MITgcmTools.read_bin
MITgcmTools.read_nctiles
MITgcmTools.MatrixInterp
— MethodMatrixInterp(in::Array{T,N},MTRX,siz) where {T,N}
Interpolate in
using MTRX
to grid of size siz
.
MITgcmTools.MetaFileRead
— MethodMetaFileRead(filIn::String)
Reads a meta file generated by MITgcm
MITgcmTools.MixedLayerDepth
— MethodMixedLayerDepth(Θ,Σ,Δ,mthd)
Compute mixed layer depth from potential temperature (Θ in °C), salinity (Σ in psu) and depth (Δ in method) according to various formulas (mthd == "BM", "Suga", "Kara"). Inputs must be dense vectors without any missing value (or NaN, etc).
D=collect(0.0:1.0:500.0); tmp=(1.0.-tanh.(5*(-1 .+ 2/D[end]*D)));
T=2.0 .+ 8.0*tmp; S=34.0 .+ 0.5*tmp;
(ρP,ρI,ρR) = SeaWaterDensity(T,S,D);
mld=MixedLayerDepth(T,S,D,"BM"); isapprox(mld,134.0)
using Plots
plot(ρP,-D,w=2,label="Potential Density",ylabel="Depth")
plot!(vec([ρP[1] ρP[end]]),-fill(mld,2),label="Mixed Layer Depth",w=2,c="black",s=:dash)
MITgcmTools.SeaWaterDensity
— FunctionSeaWaterDensity(Θ,Σ,Π,Π0)
Compute potential density (ρP), in situ density (ρI), and density referenced to PREF (Π0 in decibars) from potential temperature (Θ in °C), salinity (Σ in psu) and pressure (Π in decibars) according to the UNESCO / Jackett & McDougall 1994 equation of state.
Credits: code based on a Matlab implementation by B. Ferron Reference: https://www.jodc.go.jp/info/iocdoc/UNESCOtech/059832eb.pdf Check value: ρI = 1041.83267kg/m^3
for Θ=3°Celcius
, Σ=35psu
, Π=3000dbar
(ρP,ρI,ρR) = SeaWaterDensity(3.,35.5,3000.)
isapprox(ρI,1041.83267, rtol=1e-6)
MITgcmTools.convert2array
— Methodconvert2array(fld::MeshArray)
Convert MeshArray to Array (or vice versa otherwise)
MITgcmTools.convert2gcmfaces
— Methodconvert2gcmfaces(fld::MeshArray)
Convert mitgcm output to MeshArray (or vice versa otherwise)
MITgcmTools.findtiles
— Methodfindtiles(ni::Int,nj::Int,mygrid::gcmgrid)
findtiles(ni::Int,nj::Int,grid::String="LatLonCap",GridParentDir="../inputs/GRID_LLC90/")
Return a MeshArray
map of tile indices, mytiles["tileNo"]
, for tile size ni,nj
and extract grid variables accordingly.
MITgcmTools.parsemeta
— Methodparsemeta(metafile)
Parse out an MITgcm
metadata file and return a Dict
of fields in the file.
MITgcmTools.prep_MTRX
— Methodprep_MTRX()
Repackage interpolation matrix, mask, etc to .jld
file.
MITgcmTools.qwckplot
— Methodqwckplot(fld::MeshArray,ttl::String)
Plot input using convert2array and heatmap + add title
MITgcmTools.qwckplot
— Methodqwckplot(fld::MeshArray,clims::NTuple{2, Number})
Plot input using convert2array and heatmap w. chosen clims
MITgcmTools.qwckplot
— Methodqwckplot(fld::MeshArray)
Call qwckplot(fld::MeshArray) with date as title. Example:
!isdir("GRID_LLC90") ? error("missing files") : nothing
GridVariables=GridLoad(GridSpec("LLC90"))
qwckplot(GridVariables["Depth"])
MITgcmTools.readAvailDiagnosticsLog
— MethodreadAvailDiagnosticsLog(fname,fldname)
Get the information for a particular field from the available_diagnostics.log
file (MITgcm
output).
MITgcmTools.read_SPM
— Functionread_SPM(pth::String,fil::String="interp_precomputed.mat")
Reads pre-computed interpolation (sparse matrix) from pth*"interp_precomputed.mat"
.
MITgcmTools.read_bin
— Methodread_bin(fil::String,kt::Union{Int,Missing},kk::Union{Int,Missing},prec::DataType,mygrid::gcmgrid)
Read model output from binary file and convert to MeshArray. Other methods:
read_bin(fil::String,prec::DataType,mygrid::gcmgrid)
read_bin(fil::String,mygrid::gcmgrid)
MITgcmTools.read_nctiles
— Methodread_nctiles(fileName,fldName,mygrid)
Read model output from NCTiles file and convert to MeshArray instance.
mygrid=GridSpec("LatLonCap")
fileName="nctiles_grid/GRID"
Depth=read_nctiles(fileName,"Depth",mygrid)
hFacC=read_nctiles(fileName,"hFacC",mygrid)
hFacC=read_nctiles(fileName,"hFacC",mygrid,I=(:,:,1))