archive-gh-me/MeteoInfo
1
0
Fork 0
mirror of https://github.com/meteoinfo/MeteoInfo.git synced 2025-12-08 20:36:05 +00:00
Code Issues Projects Releases Wiki Activity

re-arrange meteolib's packages and modules

Browse Source
This commit is contained in:
wyq 2020-08-08 16:50:24 +08:00
parent cc20e1e7b9
commit 4af9c820a2
14 changed files with 483 additions and 217 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

16
.idea/inspectionProfiles/Project_Default.xml generated
View File

@ -300,22 +300,6 @@
<value />
</option>
</inspection_tool>
<inspection_tool class="PyInterpreterInspection" enabled="false" level="WARNING" enabled_by_default="false" />
<inspection_tool class="PyStatementEffectInspection" enabled="false" level="WARNING" enabled_by_default="false" />
<inspection_tool class="PyUnresolvedReferencesInspection" enabled="true" level="WARNING" enabled_by_default="true">
<option name="ignoredIdentifiers">
<list>
<option value="pylib.mipylib.numeric.core.multiarray.*" />
<option value="numbers" />
<option value="jarray" />
<option value="datetime" />
<option value="pylib.mipylib.numeric.core._dtype.*" />
<option value="pylib.mipylib.numeric.core.numeric.*" />
<option value="math" />
<option value="cmath" />
</list>
</option>
</inspection_tool>
<inspection_tool class="RawTypeCanBeGeneric" enabled="false" level="INFORMATION" enabled_by_default="false" />
<inspection_tool class="RawUseOfParameterizedType" enabled="false" level="WARNING" enabled_by_default="false" />
<inspection_tool class="ReadWriteStringCanBeUsed" enabled="false" level="WARNING" enabled_by_default="false" />

16
MeteoInfoLab/milconfig.xml
View File

@ -1,8 +1,6 @@
<?xml version="1.0" encoding="UTF-8" standalone="no"?>
<MeteoInfo File="milconfig.xml" Type="configurefile">
<Path OpenPath="D:\Working\MIScript\Jython\mis\common_math\linalg">
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\data_process"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\dataset"/>
<Path OpenPath="D:\Working\MIScript\Jython\mis\meteo\interpolation">
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\hdf"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\io\grib"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\io"/>
@ -10,25 +8,25 @@
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\common_math\interpolate"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\common_math\fitting"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\grads"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\plot_types"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\plot_types\3d"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\plot_types\3d\jogl"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\common_math"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\common_math\linalg"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\common_math"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\meteo"/>
<RecentFolder Folder="D:\Working\MIScript\Jython\mis\meteo\interpolation"/>
</Path>
<File>
<OpenedFiles>
<OpenedFile File="D:\Working\MIScript\Jython\mis\toolbox\miml\model_selection\train_test_split_2.py"/>
<OpenedFile File="D:\Working\MIScript\Jython\mis\toolbox\miml\deep_learning\mnist_1.py"/>
<OpenedFile File="D:\Working\MIScript\Jython\mis\plot_types\3d\jogl\plot_cuace_3d_particles_relief.py"/>
<OpenedFile File="D:\Working\MIScript\Jython\mis\common_math\linalg\solve_1.py"/>
<OpenedFile File="D:\Working\MIScript\Jython\mis\meteo\interpolation\log_interpolate_1d_1.py"/>
</OpenedFiles>
<RecentFiles>
<RecentFile File="D:\Working\MIScript\Jython\mis\toolbox\miml\model_selection\train_test_split_2.py"/>
<RecentFile File="D:\Working\MIScript\Jython\mis\toolbox\miml\deep_learning\mnist_1.py"/>
<RecentFile File="D:\Working\MIScript\Jython\mis\plot_types\3d\jogl\plot_cuace_3d_particles_relief.py"/>
<RecentFile File="D:\Working\MIScript\Jython\mis\common_math\linalg\solve_1.py"/>
<RecentFile File="D:\Working\MIScript\Jython\mis\meteo\interpolation\log_interpolate_1d_1.py"/>
</RecentFiles>
</File>
<Font>

BIN
MeteoInfoLab/pylib/mipylib/meteolib/__init__$py.class
View File

Binary file not shown.

14
MeteoInfoLab/pylib/mipylib/meteolib/__init__.py
View File

@ -1,5 +1,9 @@
from .meteo import *
import wrf
__all__ = []
__all__ += meteo.__all__
from .meteo import *
import wrf
from .calc import *
from .interpolate import *
__all__ = []
__all__ += meteo.__all__
__all__ += calc.__all__
__all__ += interpolate.__all__

6
MeteoInfoLab/pylib/mipylib/meteolib/calc/__init__.py Normal file
View File

@ -0,0 +1,6 @@
from .kinematics import *
from .thermo import *
__all__ = []
__all__ += kinematics.__all__
__all__ += thermo.__all__

36
MeteoInfoLab/pylib/mipylib/meteolib/calc/kinematics.py Normal file
View File

@ -0,0 +1,36 @@
"""Contains calculation of kinematic parameters (e.g. divergence or vorticity)."""
from org.meteoinfo.math.meteo import MeteoMath
import mipylib.numeric as np
from mipylib.numeric.core import NDArray, DimArray
__all__ = [
'vorticity'
]
def vorticity(u, v, x=None, y=None):
"""
Calculates the vertical component of the curl (ie, vorticity). The data should be lon/lat projection.
:param u: (*array*) U component array (2D).
:param v: (*array*) V component array (2D).
:param x: (*array*) X coordinate array (1D).
:param y: (*array*) Y coordinate array (1D).
:returns: Array of the vertical component of the curl.
"""
if x is None or y is None:
if isinstance(u, DimArray) and isinstance(v, DimArray):
x = u.dimvalue(1)
y = u.dimvalue(0)
else:
raise ValueError("Need x, y coordinates")
if isinstance(x, (list, tuple)):
x = np.array(x)
if isinstance(y, (list, tuple)):
y = np.array(y)
r = MeteoMath.hcurl(u.asarray(), v.asarray(), x.asarray(), y.asarray())
return DimArray(NDArray(r), u.dims, u.fill_value, u.proj)

130
MeteoInfoLab/pylib/mipylib/meteolib/calc/thermo.py Normal file
View File

@ -0,0 +1,130 @@
"""Contains a collection of thermodynamic calculations."""
from .. import constants
import mipylib.numeric as np
__all__ = [
'mixing_ratio','mixing_ratio_from_specific_humidity','relative_humidity_from_specific_humidity',
'saturation_mixing_ratio','saturation_vapor_pressure'
]
def saturation_vapor_pressure(temperature):
r"""Calculate the saturation water vapor (partial) pressure.
Parameters
----------
temperature : `float`
The temperature (celsius)
Returns
-------
`float`
The saturation water vapor (partial) pressure
See Also
--------
vapor_pressure, dewpoint
Notes
-----
Instead of temperature, dewpoint may be used in order to calculate
the actual (ambient) water vapor (partial) pressure.
The formula used is that from [Bolton1980]_ for T in degrees Celsius:
.. math:: 6.112 e^\frac{17.67T}{T + 243.5}
"""
# Converted from original in terms of C to use kelvin. Using raw absolute values of C in
# a formula plays havoc with units support.
return constants.sat_pressure_0c * np.exp(17.67 * temperature
/ (temperature + 243.5))
def mixing_ratio_from_specific_humidity(specific_humidity):
r"""Calculate the mixing ratio from specific humidity.
Parameters
----------
specific_humidity: `pint.Quantity`
Specific humidity of air
Returns
-------
`pint.Quantity`
Mixing ratio
Notes
-----
Formula from [Salby1996]_ pg. 118.
.. math:: w = \frac{q}{1-q}
* :math:`w` is mixing ratio
* :math:`q` is the specific humidity
See Also
--------
mixing_ratio, specific_humidity_from_mixing_ratio
"""
return specific_humidity / (1 - specific_humidity)
def mixing_ratio(part_press, tot_press):
"""
Calculates the mixing ratio of gas given its partial pressure
and the total pressure of the air.
There are no required units for the input arrays, other than that
they have the same units.
Parameters
----------
part_press : array_like
Partial pressure of the constituent gas
tot_press : array_like
Total air pressure
Returns
-------
array_like
The (mass) mixing ratio, dimensionless (e.g. Kg/Kg or g/g)
See Also
--------
vapor_pressure
"""
return constants.epsilon * part_press / (tot_press - part_press)
def saturation_mixing_ratio(tot_press, temperature):
"""
Calculates the saturation mixing ratio given total pressure
and the temperature.
The implementation uses the formula outlined in [4]
Parameters
----------
tot_press: array_like
Total atmospheric pressure (hPa)
temperature: array_like
The temperature (celsius)
Returns
-------
array_like
The saturation mixing ratio, dimensionless
References
----------
.. [4] Hobbs, Peter V. and Wallace, John M., 1977: Atmospheric Science, an Introductory
Survey. 73.
"""
return mixing_ratio(saturation_vapor_pressure(temperature), tot_press)
def relative_humidity_from_specific_humidity(specific_humidity, temperature, pressure):
r"""Calculate the relative humidity from specific humidity, temperature, and pressure.
Parameters
----------
specific_humidity: `pint.Quantity`
Specific humidity of air
temperature: `pint.Quantity`
Air temperature
pressure: `pint.Quantity`
Total atmospheric pressure
Returns
-------
`pint.Quantity`
Relative humidity
Notes
-----
Formula based on that from [Hobbs1977]_ pg. 74. and [Salby1996]_ pg. 118.
.. math:: RH = \frac{q}{(1-q)w_s}
* :math:`RH` is relative humidity as a unitless ratio
* :math:`q` is specific humidity
* :math:`w_s` is the saturation mixing ratio
See Also
--------
relative_humidity_from_mixing_ratio
"""
return (mixing_ratio_from_specific_humidity(specific_humidity)
/ saturation_mixing_ratio(pressure, temperature))

4
MeteoInfoLab/pylib/mipylib/meteolib/interpolate/__init__.py Normal file
View File

@ -0,0 +1,4 @@
from .one_dimension import *
__all__ = []
__all__ += one_dimension.__all__

75
MeteoInfoLab/pylib/mipylib/meteolib/interpolate/one_dimension.py Normal file
View File

@ -0,0 +1,75 @@
"""Interpolate data along a single axis."""
from org.meteoinfo.math import ArrayUtil
import mipylib.numeric as np
from mipylib.numeric.core import NDArray
__all__ = [
'interpolate_1d','log_interpolate_1d'
]
def interpolate_1d(x, xp, *args, **kwargs):
'''
Interpolation over a specified axis for arrays of any shape.
Parameters
----------
x : array-like
1-D array of desired interpolated values.
xp : array-like
The x-coordinates of the data points.
args : array-like
The data to be interpolated. Can be multiple arguments, all must be the same shape as
xp.
axis : int, optional
The axis to interpolate over. Defaults to 0.
Returns
-------
array-like
Interpolated values for each point with coordinates sorted in ascending order.
'''
axis = kwargs.pop('axis', 0)
if isinstance(x, (list, tuple)):
x = np.array(x)
if isinstance(x, NDArray):
x = x._array
vars = args
ret = []
for a in vars:
r = ArrayUtil.interpolate_1d(x, xp._array, a._array, axis)
ret.append(NDArray(r))
if len(ret) == 1:
return ret[0]
else:
return ret
def log_interpolate_1d(x, xp, *args, **kwargs):
'''
Interpolation on a logarithmic x-scale for interpolation values in pressure coordintates.
Parameters
----------
x : array-like
1-D array of desired interpolated values.
xp : array-like
The x-coordinates of the data points.
args : array-like
The data to be interpolated. Can be multiple arguments, all must be the same shape as
xp.
axis : int, optional
The axis to interpolate over. Defaults to 0.
Returns
-------
array-like
Interpolated values for each point with coordinates sorted in ascending order.
'''
# Log x and xp
log_x = np.log(x)
log_xp = np.log(xp)
return interpolate_1d(log_x, log_xp, *args, **kwargs)

BIN
MeteoInfoLab/pylib/mipylib/meteolib/meteo$py.class
View File

Binary file not shown.

214
MeteoInfoLab/pylib/mipylib/meteolib/meteo.py
View File

@ -16,11 +16,10 @@ import constants as constants
__all__ = [
'cumsimp','dewpoint','dewpoint2rh','dewpoint_rh','dry_lapse','ds2uv','equivalent_potential_temperature',
'exner_function','flowfun','h2p',
'interpolate_1d','log_interpolate_1d','mixing_ratio','mixing_ratio_from_specific_humidity',
'moist_lapse','p2h','potential_temperature','qair2rh','rh2dewpoint','relative_humidity_from_specific_humidity',
'saturation_mixing_ratio','saturation_vapor_pressure','sigma_to_pressure','tc2tf',
'moist_lapse','p2h','potential_temperature','qair2rh','rh2dewpoint',
'sigma_to_pressure','tc2tf',
'temperature_from_potential_temperature','tf2tc','uv2ds','pressure_to_height_std',
'height_to_pressure_std','eof','vapor_pressure','varimax'
'height_to_pressure_std','eof','vapor_pressure','varimax','virtual_temperature'
]
def uv2ds(u, v):
@ -239,56 +238,6 @@ def dewpoint2rh(dewpoint, temp):
else:
return MeteoMath.dewpoint2rh(temp, dewpoint)
def mixing_ratio_from_specific_humidity(specific_humidity):
r"""Calculate the mixing ratio from specific humidity.
Parameters
----------
specific_humidity: `pint.Quantity`
Specific humidity of air
Returns
-------
`pint.Quantity`
Mixing ratio
Notes
-----
Formula from [Salby1996]_ pg. 118.
.. math:: w = \frac{q}{1-q}
* :math:`w` is mixing ratio
* :math:`q` is the specific humidity
See Also
--------
mixing_ratio, specific_humidity_from_mixing_ratio
"""
return specific_humidity / (1 - specific_humidity)
def relative_humidity_from_specific_humidity(specific_humidity, temperature, pressure):
r"""Calculate the relative humidity from specific humidity, temperature, and pressure.
Parameters
----------
specific_humidity: `pint.Quantity`
Specific humidity of air
temperature: `pint.Quantity`
Air temperature
pressure: `pint.Quantity`
Total atmospheric pressure
Returns
-------
`pint.Quantity`
Relative humidity
Notes
-----
Formula based on that from [Hobbs1977]_ pg. 74. and [Salby1996]_ pg. 118.
.. math:: RH = \frac{q}{(1-q)w_s}
* :math:`RH` is relative humidity as a unitless ratio
* :math:`q` is specific humidity
* :math:`w_s` is the saturation mixing ratio
See Also
--------
relative_humidity_from_mixing_ratio
"""
return (mixing_ratio_from_specific_humidity(specific_humidity)
/ saturation_mixing_ratio(pressure, temperature))
def rh2dewpoint(rh, temp):
"""
Calculate dewpoint from relative humidity and temperature
@ -452,53 +401,7 @@ def moist_lapse(pressure, temperature):
(constants.Cp_d + (constants.Lv * constants.Lv * rs * constants.epsilon / (constants.Rd * t * t)))).to('kelvin')
return frac / p
return dt
def mixing_ratio(part_press, tot_press):
"""
Calculates the mixing ratio of gas given its partial pressure
and the total pressure of the air.
There are no required units for the input arrays, other than that
they have the same units.
Parameters
----------
part_press : array_like
Partial pressure of the constituent gas
tot_press : array_like
Total air pressure
Returns
-------
array_like
The (mass) mixing ratio, dimensionless (e.g. Kg/Kg or g/g)
See Also
--------
vapor_pressure
"""
return constants.epsilon * part_press / (tot_press - part_press)
def saturation_mixing_ratio(tot_press, temperature):
"""
Calculates the saturation mixing ratio given total pressure
and the temperature.
The implementation uses the formula outlined in [4]
Parameters
----------
tot_press: array_like
Total atmospheric pressure (hPa)
temperature: array_like
The temperature (celsius)
Returns
-------
array_like
The saturation mixing ratio, dimensionless
References
----------
.. [4] Hobbs, Peter V. and Wallace, John M., 1977: Atmospheric Science, an Introductory
Survey. 73.
"""
return mixing_ratio(saturation_vapor_pressure(temperature), tot_press)
def vapor_pressure(pressure, mixing):
r"""Calculate water vapor (partial) pressure.
Given total `pressure` and water vapor `mixing` ratio, calculates the
@ -524,32 +427,7 @@ def vapor_pressure(pressure, mixing):
saturation_vapor_pressure, dewpoint
"""
return pressure * mixing / (constants.epsilon + mixing)
def saturation_vapor_pressure(temperature):
r"""Calculate the saturation water vapor (partial) pressure.
Parameters
----------
temperature : `float`
The temperature (celsius)
Returns
-------
`float`
The saturation water vapor (partial) pressure
See Also
--------
vapor_pressure, dewpoint
Notes
-----
Instead of temperature, dewpoint may be used in order to calculate
the actual (ambient) water vapor (partial) pressure.
The formula used is that from [Bolton1980]_ for T in degrees Celsius:
.. math:: 6.112 e^\frac{17.67T}{T + 243.5}
"""
# Converted from original in terms of C to use kelvin. Using raw absolute values of C in
# a formula plays havoc with units support.
return constants.sat_pressure_0c * np.exp(17.67 * temperature
/ (temperature + 243.5))
def exner_function(pressure, reference_pressure=constants.P0):
r"""Calculate the Exner function.
.. math:: \Pi = \left( \frac{p}{p_0} \right)^\kappa
@ -647,72 +525,36 @@ def temperature_from_potential_temperature(pressure, theta):
>>> T = temperature_from_potential_temperature(p,theta)
"""
return theta * exner_function(pressure)
def interpolate_1d(x, xp, *args, **kwargs):
'''
Interpolation over a specified axis for arrays of any shape.
def virtual_temperature(temperature, mixing, molecular_weight_ratio=constants.epsilon):
r"""Calculate virtual temperature.
This calculation must be given an air parcel's temperature and mixing ratio.
The implementation uses the formula outlined in [Hobbs2006]_ pg.80.
Parameters
----------
x : array-like
1-D array of desired interpolated values.
xp : array-like
The x-coordinates of the data points.
args : array-like
The data to be interpolated. Can be multiple arguments, all must be the same shape as
xp.
axis : int, optional
The axis to interpolate over. Defaults to 0.
temperature: `array`
air temperature
mixing : `array`
dimensionless mass mixing ratio
molecular_weight_ratio : float, optional
The ratio of the molecular weight of the constituent gas to that assumed
for air. Defaults to the ratio for water vapor to dry air.
(:math:`\epsilon\approx0.622`).
Returns
-------
array-like
Interpolated values for each point with coordinates sorted in ascending order.
'''
axis = kwargs.pop('axis', 0)
if isinstance(x, (list, tuple)):
x = np.array(x)
if isinstance(x, NDArray):
x = x._array
vars = args
ret = []
for a in vars:
r = ArrayUtil.interpolate_1d(x, xp._array, a._array, axis)
ret.append(NDArray(r))
if len(ret) == 1:
return ret[0]
else:
return ret
def log_interpolate_1d(x, xp, *args, **kwargs):
'''
Interpolation on a logarithmic x-scale for interpolation values in pressure coordintates.
Parameters
----------
x : array-like
1-D array of desired interpolated values.
xp : array-like
The x-coordinates of the data points.
args : array-like
The data to be interpolated. Can be multiple arguments, all must be the same shape as
xp.
axis : int, optional
The axis to interpolate over. Defaults to 0.
`array`
The corresponding virtual temperature of the parcel
Returns
-------
array-like
Interpolated values for each point with coordinates sorted in ascending order.
'''
# Log x and xp
log_x = np.log(x)
log_xp = np.log(xp)
return interpolate_1d(log_x, log_xp, *args, **kwargs)
Notes
-----
.. math:: T_v = T \frac{\text{w} + \epsilon}{\epsilon\,(1 + \text{w})}
"""
return temperature * ((mixing + molecular_weight_ratio)
/ (molecular_weight_ratio * (1 + mixing)))
def cumsimp(y):
"""

187
MeteoInfoLib/src/main/java/org/meteoinfo/math/meteo/MeteoMath.java
View File

@ -6,11 +6,15 @@
package org.meteoinfo.math.meteo;
import org.meteoinfo.math.ArrayMath;
import org.meteoinfo.math.ArrayUtil;
import org.meteoinfo.ndarray.Array;
import org.meteoinfo.ndarray.Index;
import org.meteoinfo.ndarray.DataType;
import org.meteoinfo.ndarray.IndexIterator;
import java.util.Arrays;
import java.util.List;
/**
*
* @author yaqiang
@ -915,4 +919,187 @@ public class MeteoMath {
return new double[]{u, v};
}
/**
* Performs a centered difference operation on a grid data along one
* dimension direction
*
* @param data The grid data
* @param dimIdx Direction dimension index
* @return Result grid data
*/
public static Array cdiff(Array data, int dimIdx) {
Array r = Array.factory(DataType.DOUBLE, data.getShape());
Index index = data.getIndex();
Index indexr = r.getIndex();
int[] shape = data.getShape();
int[] current, cc;
double a, b;
for (int i = 0; i < r.getSize(); i++) {
current = indexr.getCurrentCounter();
if (current[dimIdx] == 0 || current[dimIdx] == shape[dimIdx] - 1) {
r.setDouble(indexr, Double.NaN);
} else {
cc = Arrays.copyOf(current, current.length);
cc[dimIdx] = cc[dimIdx] - 1;
index.set(cc);
a = data.getDouble(index);
cc[dimIdx] = cc[dimIdx] + 2;
index.set(cc);
b = data.getDouble(index);
if (Double.isNaN(a) || Double.isNaN(b)) {
r.setDouble(indexr, Double.NaN);
} else {
r.setDouble(indexr, a - b);
}
}
indexr.incr();
}
return r;
}
/**
* Performs a centered difference operation on a grid data in the x or y
* direction
*
* @param data The grid data
* @param isX If is x direction
* @return Result grid data
*/
public static Array cdiff_bak(Array data, boolean isX) {
if (data.getRank() == 2) {
int xnum = data.getShape()[1];
int ynum = data.getShape()[0];
Array r = Array.factory(DataType.DOUBLE, data.getShape());
for (int i = 0; i < ynum; i++) {
for (int j = 0; j < xnum; j++) {
if (i == 0 || i == ynum - 1 || j == 0 || j == xnum - 1) {
r.setDouble(i * xnum + j, Double.NaN);
} else {
double a, b;
if (isX) {
a = data.getDouble(i * xnum + j + 1);
b = data.getDouble(i * xnum + j - 1);
} else {
a = data.getDouble((i + 1) * xnum + j);
b = data.getDouble((i - 1) * xnum + j);
}
if (Double.isNaN(a) || Double.isNaN(b)) {
r.setDouble(i * xnum + j, Double.NaN);
} else {
r.setDouble(i * xnum + j, a - b);
}
}
}
}
return r;
} else if (data.getRank() == 1) {
int n = data.getShape()[0];
Array r = Array.factory(DataType.DOUBLE, data.getShape());
for (int i = 0; i < n; i++) {
if (i == 0 || i == n - 1) {
r.setDouble(i, Double.NaN);
} else {
double a, b;
a = data.getDouble(i + 1);
b = data.getDouble(i - 1);
if (Double.isNaN(a) || Double.isNaN(b)) {
r.setDouble(i, Double.NaN);
} else {
r.setDouble(i, a - b);
}
}
}
return r;
} else {
System.out.println("Data dimension number must be 1 or 2!");
return null;
}
}
/**
* Calculates the vertical component of the curl (ie, vorticity)
*
* @param uData U component
* @param vData V component
* @param xx X dimension value
* @param yy Y dimension value
* @return Curl
*/
public static Array hcurl(Array uData, Array vData, Array xx, Array yy) {
int rank = uData.getRank();
int[] shape = uData.getShape();
Array[] llData = ArrayUtil.meshgrid(xx, yy);
Array lonData = llData[0];
Array latData = llData[1];
Array dv = cdiff(vData, rank - 1);
Array dx = ArrayMath.mul(cdiff(lonData, rank - 1), Math.PI / 180);
Array du = cdiff(ArrayMath.mul(uData, ArrayMath.cos(ArrayMath.mul(latData, Math.PI / 180))), rank - 2);
Array dy = ArrayMath.mul(cdiff(latData, rank - 2), Math.PI / 180);
Array gData = ArrayMath.div(ArrayMath.sub(ArrayMath.div(dv, dx), ArrayMath.div(du, dy)), ArrayMath.mul(ArrayMath.cos(ArrayMath.mul(latData, Math.PI / 180)), 6.37e6));
return gData;
}
/**
* Calculates the horizontal divergence using finite differencing
*
* @param uData U component
* @param vData V component
* @param xx X dimension value
* @param yy Y dimension value
* @return Divergence
*/
public static Array hdivg(Array uData, Array vData, List<Number> xx, List<Number> yy) {
int rank = uData.getRank();
int[] shape = uData.getShape();
Array lonData = Array.factory(DataType.DOUBLE, shape);
Array latData = Array.factory(DataType.DOUBLE, shape);
Index index = lonData.getIndex();
int[] current;
for (int i = 0; i < lonData.getSize(); i++) {
current = index.getCurrentCounter();
lonData.setDouble(index, xx.get(current[rank - 1]).doubleValue());
latData.setDouble(index, yy.get(current[rank - 2]).doubleValue());
index.incr();
}
Array du = cdiff(uData, rank - 1);
Array dx = ArrayMath.mul(cdiff(lonData, rank - 1), Math.PI / 180);
Array dv = cdiff(ArrayMath.mul(vData, ArrayMath.cos(ArrayMath.mul(latData, Math.PI / 180))), rank - 2);
Array dy = ArrayMath.mul(cdiff(latData, rank - 2), Math.PI / 180);
Array gData = ArrayMath.div(ArrayMath.add(ArrayMath.div(du, dx), ArrayMath.div(dv, dy)), ArrayMath.mul(ArrayMath.cos(ArrayMath.mul(latData, Math.PI / 180)), 6.37e6));
return gData;
}
/**
* Take magnitude value from U/V grid data
*
* @param uData U grid data
* @param vData V grid data
* @return Magnitude grid data
*/
public static Array magnitude(Array uData, Array vData) {
int[] shape = uData.getShape();
Array r = Array.factory(DataType.DOUBLE, shape);
IndexIterator iterU = uData.getIndexIterator();
IndexIterator iterV = vData.getIndexIterator();
IndexIterator iterR = r.getIndexIterator();
double u, v;
while (iterU.hasNext()) {
u = iterU.getDoubleNext();
v = iterV.getDoubleNext();
if (Double.isNaN(u) || Double.isNaN(v)) {
iterR.setDoubleNext(Double.NaN);
} else {
iterR.setDoubleNext(Math.sqrt(Math.pow(u, 2) + Math.pow(v, 2)));
}
}
return r;
}
}

BIN
MeteoInfoLib/target/classes/org/meteoinfo/math/meteo/MeteoMath.class
View File

Binary file not shown.

2
MeteoInfoMap/config.xml
View File

@ -1,6 +1,6 @@
<?xml version="1.0" encoding="UTF-8" standalone="no"?>
<MeteoInfo File="config.xml" Type="configurefile">
<Path OpenPath="D:\Temp\micaps"/>
<Path OpenPath="D:\Temp\test"/>
<Font>
<TextFont FontName="Arial" FontSize="14"/>
<LegendFont FontName="宋体" FontSize="12"/>