libGM76.pl
changeset 39 56bdfe65a697
parent 38 15c603bc4f70
equal deleted inserted replaced
38:15c603bc4f70 39:56bdfe65a697
     1 #======================================================================
     1 #======================================================================
     2 #                    L I B G M 7 6 . P L 
     2 #                    L I B G M 7 6 . P L 
     3 #                    doc: Sun Feb 20 14:43:47 2011
     3 #                    doc: Sun Feb 20 14:43:47 2011
     4 #                    dlm: Sat Feb  2 02:33:11 2019
     4 #                    dlm: Sat Apr  6 19:33:05 2019
     5 #                    (c) 2011 A.M. Thurnherr
     5 #                    (c) 2011 A.M. Thurnherr
     6 #                    uE-Info: 144 43 NIL 0 0 70 2 2 4 NIL ofnI
     6 #                    uE-Info: 34 48 NIL 0 0 70 2 2 4 NIL ofnI
     7 #======================================================================
     7 #======================================================================
     8 
     8 
     9 # HISTORY:
     9 # HISTORY:
    10 #	Feb 20, 2011: - created
    10 #	Feb 20, 2011: - created
    11 #	Feb 28, 2011: - cosmetics
    11 #	Feb 28, 2011: - cosmetics
    25 #				  - replaced old code based on Gregg + Kunze formalism with
    25 #				  - replaced old code based on Gregg + Kunze formalism with
    26 #					expressions from Thurnherr et al. (GRL 2015)
    26 #					expressions from Thurnherr et al. (GRL 2015)
    27 #				  - added GM_strain
    27 #				  - added GM_strain
    28 #				  - BUG: Sw usage message had wrong parameter order
    28 #				  - BUG: Sw usage message had wrong parameter order
    29 #				  - renamed Sw => GM_VKE; Su_z => GM_shear
    29 #				  - renamed Sw => GM_VKE; Su_z => GM_shear
    30 
    30 #	Mar 31, 2019: - updated doc for shear and strain spectra
       
    31 #	Apr  1, 2019: - added GM_vdiv
       
    32 #	Apr  5, 2019: - BUG: GM_VKE was erroneous (argument shifting was wrong)
       
    33 #				  - adapted to improved antsFunUsage()
       
    34 #	Apr  6, 2019: - fiddled during debugging of fs_finestructure issue
    31 
    35 
    32 require "$ANTS/libEOS83.pl";
    36 require "$ANTS/libEOS83.pl";
    33 
    37 
    34 my($pi) = 3.14159265358979;
    38 my($pi) = 3.14159265358979;
    35 
    39 
    98 }
   102 }
    99 
   103 
   100 sub GM_VKE(@)
   104 sub GM_VKE(@)
   101 {
   105 {
   102 	my($omega,$m,$lat,$b,$jstar,$N) =
   106 	my($omega,$m,$lat,$b,$jstar,$N) =
   103 		&antsFunUsage(-5,'fff','[frequency[1/s]] <vertical wavenumber[rad/m]> <lat[deg]> <b[m]> <j*> <N[rad/s]>',@_);
   107 		&antsFunUsage(-5,'fffff','GM_VKE([frequency[1/s]] <vertical wavenumber[rad/m]> <lat[deg]> <b[m]> <j*> <N[rad/s]>)',@_);
   104 
   108 
   105 	if (defined($N)) {									# Sw(omega,m)
   109 	if (defined($N)) {									# Sw(omega,m)
   106 		local($f) = abs(&f($lat));
   110 		local($f) = abs(&f($lat));
   107 		$omega += $PRACTICALLY_ZERO if ($omega < $f);
   111 		$omega += $PRACTICALLY_ZERO if ($omega < $f);
   108 		$omega -= $PRACTICALLY_ZERO if ($omega > $N);
   112 		$omega -= $PRACTICALLY_ZERO if ($omega > $N);
   109 		return nan if ($omega < $f || $omega > $N);
   113 		return nan if ($omega < $f || $omega > $N);
   110 		my($mstar) = &m($jstar,$N,$omega);
   114 		my($mstar) = &m($jstar,$N,$omega);
   111 		return $GM_E0 * $b * 2 * $f**2/$omega**2/B($omega) * $jstar / ($m+$mstar)**2;
   115 		return $GM_E0 * $b * 2 * $f**2/$omega**2/B($omega) * $jstar / ($m+$mstar)**2;
   112 	} else {											# Sw(m), i.e. integrated over all omega; as in Thurnherr et al., GRL 2015
   116 	} else {											# Sw(m), i.e. integrated over all omega; as in Thurnherr et al., GRL 2015
   113 		$N = $lat;										# shift arguments to account for missing omega
   117 		$N = $jstar;  $jstar = $b;						# shift arguments to account for missing omega
   114 		$lat = $m;
   118 		$b = $lat; $lat = $m;
   115 		local($f) = abs(&f($lat));
   119 		local($f) = abs(&f($lat));
   116 		$m = $omega;
   120 		$m = $omega;
   117 		undef($omega);
   121 		undef($omega);
   118 		my($mstar) = &m($jstar,$N);
   122 		my($mstar) = &m($jstar,$N);
   119 		return $pi * $GM_E0 * $b * $N * $f * $jstar / ($m+$mstar)**2;
   123 		return $pi * $GM_E0 * $b * $N * $f * $jstar / ($m+$mstar)**2;
   120 	}
   124 	}
   121 }
   125 }
   122 
   126 
   123 #----------------------------------------------------------------------
   127 #----------------------------------------------------------------------
       
   128 # Vertical Divergence (normalized by N)
       
   129 #	- implemented to investigate shear-to-vertical divergence ratio
       
   130 #	- S[w_z/N] = S[w_z] / N^2 = S[w] * m^2 / N^2
       
   131 #	- GM shear-to-vd variance ratio = 3N/2f
       
   132 #	- GM strain-to-vd variance ratio = N/2f
       
   133 #----------------------------------------------------------------------
       
   134 
       
   135 sub GM_vdiv(@)
       
   136 {
       
   137 	my($m,$lat,$b,$jstar,$N) =
       
   138 		&antsFunUsage(5,'fffff','GM_vdiv(<vertical wavenumber[rad/m]> <lat[deg]> <b[m]> <j*> <N[rad/s]>)',@_);
       
   139 	return GM_VKE($m,$lat,$b,$jstar,$N) * $m**2 / $N**2;
       
   140 }
       
   141 
       
   142 #----------------------------------------------------------------------
   124 # Shear and Strain m-spectra (i.e. integrated over f)
   143 # Shear and Strain m-spectra (i.e. integrated over f)
       
   144 #	- shear is buoyancy-frequency normalized
   125 #	- spectral density
   145 #	- spectral density
   126 #	- from Thurnherr et al. (GRL 2015)
   146 #	- from Thurnherr et al. (GRL 2015), which is from email info
       
   147 #	  Eric Kunze, 04/22/2015: First, the standard GHP parameterization that I
       
   148 #	  implement uses either the N-normalized shear or strain spectra. The GM
       
   149 #	  versions of these are
       
   150 #		S[V_z/N](m) = (3*PI/2)*E0*b*jstar*m^2/(m+mstar)^2 and
       
   151 #		S[Z_z](m) = (PI/2)*E0*b*jstar*m^2/(m+mstar)^2
       
   152 #	  though in practice mstar is not particularly important because the
       
   153 #	  variance of both these quantities is dominated by high m.
   127 #----------------------------------------------------------------------
   154 #----------------------------------------------------------------------
   128 
   155 
   129 sub GM_shear(@)
   156 sub GM_shear(@)
   130 {
   157 {
   131 	my($m,$lat,$b,$jstar,$N) =
   158 	my($m,$lat,$b,$jstar,$N) =
   132 		&antsFunUsage(5,'fffff','<vertical wavenumber[rad/m]> <lat[deg]> <b[m]> <j*> <N[rad/s]>',@_);
   159 		&antsFunUsage(5,'fffff','GM_shear(<vertical wavenumber[rad/m]> <lat[deg]> <b[m]> <j*> <N[rad/s]>)',@_);
   133 	local($f) = abs(&f($lat));
   160 	local($f) = abs(&f($lat));
   134 	my($mstar) = &m($jstar,$N);
   161 	my($mstar) = &m($jstar,$N);
   135 	return 3 * $pi/2 * $GM_E0 * $b * $jstar * $m**2 / ($m+$mstar)**2;
   162 	return 3 * $pi/2 * $GM_E0 * $b * $jstar * $m**2 / ($m+$mstar)**2;
   136 }
   163 }
   137 
   164 
   138 sub GM_strain(@)
   165 sub GM_strain(@)
   139 {
   166 {
   140 	my($m,$lat,$b,$jstar,$N) =
   167 	my($m,$lat,$b,$jstar,$N) =
   141 		&antsFunUsage(5,'fffff','<vertical wavenumber[rad/m]> <lat[deg]> <b[m]> <j*> <N[rad/s]>',@_);
   168 		&antsFunUsage(5,'fffff','GM_strain(<vertical wavenumber[rad/m]> <lat[deg]> <b[m]> <j*> <N[rad/s]>)',@_);
   142 	local($f) = abs(&f($lat));
   169 	local($f) = abs(&f($lat));
   143 	my($mstar) = &m($jstar,$N);
   170 	my($mstar) = &m($jstar,$N);
   144 	return $pi/2 * $GM_E0 * $b * $jstar * $m**2 / ($m+$mstar)**2;
   171 	return $pi/2 * $GM_E0 * $b * $jstar * $m**2 / ($m+$mstar)**2;
   145 }
   172 }
   146 
   173