NumRu::GAnalysis::QG_sphere

QG on sphere with non-divergent but inaccurate geostrophic wind

Public Class Methods

div_waf(*args) click to toggle source

divergence of wave activity flux (redirected to ((<div_waf>)))

# File ../../lib/numru/ganalysis/qg.rb, line 537
def self.div_waf(*args)
  super(*args)    # defined in QG_common
end

Public Instance Methods

cos_phi(gphys) click to toggle source

cosine of latitude

# File ../../lib/numru/ganalysis/qg.rb, line 495
def cos_phi(gphys)
  lam, phi, = Planet::get_lambda_phi(gphys)
  phi.cos
end
div_h(fx, fy) click to toggle source

horizontal divergence (spherical)

# File ../../lib/numru/ganalysis/qg.rb, line 528
def div_h(fx, fy)
  Planet::div_s(fx, fy)
end
gph2psi(gph, gpref) click to toggle source

geopotential height -> QG stream function

# File ../../lib/numru/ganalysis/qg.rb, line 460
def gph2psi(gph, gpref)
  gpd = gph * Met::g - gpref
  f = f_mask0(gph)
  psi = gpd / f
  psi.name = "psi"
  psi.long_name = "QG stream function"
  psi
end
gph2psi_gpref(gph) click to toggle source

geopotential height -> QG stream function and the reference geopotential

# File ../../lib/numru/ganalysis/qg.rb, line 450
def gph2psi_gpref(gph)
  gpd, gpref = gph2gpd_gpref(gph)
  f = f_mask0(gph)
  psi = gpd / f
  psi.name = "psi"
  psi.long_name = "QG stream function"
  [psi, gpref]
end
gph2q(gph) click to toggle source

geopotential height to quasi-geostrophic potential vorticity (QGPV)

# File ../../lib/numru/ganalysis/qg.rb, line 430
def gph2q(gph)
  psi, gpref = gph2psi_gpref(gph)
  n2 = gpref2n2(gpref)
  psi2q(psi, n2)
end
gph2qb(gph) click to toggle source

same as gph2q, but the QGPV is extended to reflect the lowest-level temperature anomalies

# File ../../lib/numru/ganalysis/qg.rb, line 437
def gph2qb(gph)
  psi, gpref = gph2psi_gpref(gph)
  n2 = gpref2n2(gpref)
  psi2qb(psi, n2)
end
gph2ug_vg(gph) click to toggle source

geopotential height -> geostrophic winds

# File ../../lib/numru/ganalysis/qg.rb, line 444
def gph2ug_vg(gph)
  psi, gpref = gph2psi_gpref(gph)
  psi2ug_vg(psi)
end
grad_h(gphys) click to toggle source

horizontal gradient (spherical)

# File ../../lib/numru/ganalysis/qg.rb, line 523
def grad_h(gphys)
  Planet::grad_s(gphys)
end
psi2q(psi, n2, perturbation=false) click to toggle source

QG stream function -> QGPV

# File ../../lib/numru/ganalysis/qg.rb, line 470
def psi2q(psi, n2, perturbation=false)
  ug, vg = psi2ug_vg(psi)

  if !perturbation
    avor = Planet::absvor_s(ug,vg)
    avor.name = "qgavor"
    avor.long_name = "QG abs vor"
  else
    vor = Planet::vor_s(ug,vg)
    vor.name = "qgvor"
    vor.long_name = "QG vorticity"
    avor = vor
  end

  f = f_mask0(psi)
  qzz = gpd2qzz(psi, n2) * (f*f)

  q = avor + qzz
  q.name = "q"
  q.long_name = "QG PV"

  [q, avor, qzz, ug, vg]
end
psi2qb(psi, n2, perturbation=false) click to toggle source

same as psi2q, but the QGPV is extended to reflect the lowest-level temperature anomalies

# File ../../lib/numru/ganalysis/qg.rb, line 502
def psi2qb(psi, n2, perturbation=false)
  psie = extend_bottom(psi, nil)
  n2e = extend_bottom(n2, nil)
  results = psi2q(psie, n2e, perturbation)
  results.collect{|z| cut_bottom(z)}
end
psi2ug_vg(psi) click to toggle source

QG stream function -> geostrophic winds

# File ../../lib/numru/ganalysis/qg.rb, line 510
def psi2ug_vg(psi)
  f = f_mask0(psi)
  gpx, gpy = Planet::grad_s(psi)
  vg = gpx
  ug = -gpy
  ug.name = "ug"
  vg.name = "vg"
  ug.long_name = "ug"
  vg.long_name = "vg"
  [ug, vg]
end
waf_plumb1986_B1(psi, n2) click to toggle source

Flux of the pseudo-momentum in y direction by Plumb (1986). Specifically, B_1j in Eq.(2.9), but without the factor p. This flux is relative to the mean flow. Averaged over time (if the data is 4D).

# File ../../lib/numru/ganalysis/qg.rb, line 573
def waf_plumb1986_B1(psi, n2)
  psi_x, psi_y = Planet::grad_s(psi)
  psi_z = LogP.pcdata_dz( psi )
  f2 = f_mask0(psi) ** 2
  cosphi = cos_phi(psi)
  fx = -psi_x * psi_y * cosphi
  fy = (psi_x**2 - psi_y**2 + psi_z**2 * f2 / n2) * cosphi / 2.0
  fz = -psi_y * psi_z * (cosphi * f2) / n2
  fx.name = "B_11"
  fy.name = "B_12"
  fz.name = "B_13"
  fx.long_name = "x-comp of Px flux Plumb86"
  fy.long_name = "y-comp of Px flux Plumb86"
  fz.long_name = "z-comp of Px flux Plumb86"
  if psi.rank >= 4
    fx = fx.mean(-1)
    fy = fy.mean(-1)
    fz = fz.mean(-1)
  end
  [fx, fy, fz]
end
waf_plumb1986_B2(psi, n2) click to toggle source

Flux of the pseudo-momentum in x direction by Plumb (1986). Specifically, B_2j in Eq.(2.9), but without the factor p. This flux is relative to the mean flow. Averaged over time (if the data is 4D).

# File ../../lib/numru/ganalysis/qg.rb, line 546
def waf_plumb1986_B2(psi, n2)
  psi_x, psi_y = Planet::grad_s(psi)
  psi_z = LogP.pcdata_dz( psi )
  f2 = f_mask0(psi) ** 2
  cosphi = cos_phi(psi)
  fx = (psi_x**2 - psi_y**2 - psi_z**2 * f2 / n2) * cosphi / 2.0
  fy = psi_x * psi_y * cosphi 
  fz = psi_x * psi_z * (cosphi * f2) / n2
  fx.name = "B21"
  fy.name = "B22"
  fz.name = "B23"
  fx.long_name = "WAF B2x"
  fy.long_name = "WAF B2y"
  fz.long_name = "WAF B2z"
  if psi.rank >= 4
    fx = fx.mean(-1)
    fy = fy.mean(-1)
    fz = fz.mean(-1)
  end
  [fx, fy, fz]
end

[Validate]

Generated with the Darkfish Rdoc Generator 2.