|
1 | | -from parcels.loggers import logger |
| 1 | +from parcels.tools.loggers import logger |
| 2 | +from parcels.tools.converters import unitconverters_map, UnitConverter, Geographic, GeographicPolar |
| 3 | +from parcels.tools.error import FieldSamplingError, TimeExtrapolationError |
2 | 4 | from collections import Iterable |
3 | 5 | from py import path |
4 | 6 | import numpy as np |
5 | 7 | from ctypes import Structure, c_int, c_float, POINTER, pointer |
6 | 8 | import xarray as xr |
7 | | -from math import cos, pi |
| 9 | +from math import cos |
8 | 10 | import datetime |
9 | 11 | import math |
10 | 12 | from .grid import (RectilinearZGrid, RectilinearSGrid, CurvilinearZGrid, |
11 | 13 | CurvilinearSGrid, CGrid, GridCode) |
12 | 14 |
|
13 | 15 |
|
14 | | -__all__ = ['Field', 'VectorField', 'SummedField', 'SummedVectorField', |
15 | | - 'Geographic', 'GeographicPolar', 'GeographicSquare', 'GeographicPolarSquare'] |
16 | | - |
17 | | - |
18 | | -class FieldSamplingError(RuntimeError): |
19 | | - """Utility error class to propagate erroneous field sampling""" |
20 | | - |
21 | | - def __init__(self, x, y, z, field=None): |
22 | | - self.field = field |
23 | | - self.x = x |
24 | | - self.y = y |
25 | | - self.z = z |
26 | | - message = "%s sampled at (%f, %f, %f)" % ( |
27 | | - field.name if field else "Field", self.x, self.y, self.z |
28 | | - ) |
29 | | - super(FieldSamplingError, self).__init__(message) |
30 | | - |
31 | | - |
32 | | -class TimeExtrapolationError(RuntimeError): |
33 | | - """Utility error class to propagate erroneous time extrapolation sampling""" |
34 | | - |
35 | | - def __init__(self, time, field=None): |
36 | | - if field is not None and field.grid.time_origin and time is not None: |
37 | | - time = field.grid.time_origin + np.timedelta64(int(time), 's') |
38 | | - message = "%s sampled outside time domain at time %s." % ( |
39 | | - field.name if field else "Field", time) |
40 | | - message += " Try setting allow_time_extrapolation to True" |
41 | | - super(TimeExtrapolationError, self).__init__(message) |
42 | | - |
43 | | - |
44 | | -class UnitConverter(object): |
45 | | - """ Interface class for spatial unit conversion during field sampling |
46 | | - that performs no conversion. |
47 | | - """ |
48 | | - source_unit = None |
49 | | - target_unit = None |
50 | | - |
51 | | - def to_target(self, value, x, y, z): |
52 | | - return value |
53 | | - |
54 | | - def ccode_to_target(self, x, y, z): |
55 | | - return "1.0" |
56 | | - |
57 | | - def to_source(self, value, x, y, z): |
58 | | - return value |
59 | | - |
60 | | - def ccode_to_source(self, x, y, z): |
61 | | - return "1.0" |
62 | | - |
63 | | - |
64 | | -class Geographic(UnitConverter): |
65 | | - """ Unit converter from geometric to geographic coordinates (m to degree) """ |
66 | | - source_unit = 'm' |
67 | | - target_unit = 'degree' |
68 | | - |
69 | | - def to_target(self, value, x, y, z): |
70 | | - return value / 1000. / 1.852 / 60. |
71 | | - |
72 | | - def to_source(self, value, x, y, z): |
73 | | - return value * 1000. * 1.852 * 60. |
74 | | - |
75 | | - def ccode_to_target(self, x, y, z): |
76 | | - return "(1.0 / (1000.0 * 1.852 * 60.0))" |
77 | | - |
78 | | - def ccode_to_source(self, x, y, z): |
79 | | - return "(1000.0 * 1.852 * 60.0)" |
80 | | - |
81 | | - |
82 | | -class GeographicPolar(UnitConverter): |
83 | | - """ Unit converter from geometric to geographic coordinates (m to degree) |
84 | | - with a correction to account for narrower grid cells closer to the poles. |
85 | | - """ |
86 | | - source_unit = 'm' |
87 | | - target_unit = 'degree' |
88 | | - |
89 | | - def to_target(self, value, x, y, z): |
90 | | - return value / 1000. / 1.852 / 60. / cos(y * pi / 180) |
91 | | - |
92 | | - def to_source(self, value, x, y, z): |
93 | | - return value * 1000. * 1.852 * 60. * cos(y * pi / 180) |
94 | | - |
95 | | - def ccode_to_target(self, x, y, z): |
96 | | - return "(1.0 / (1000. * 1.852 * 60. * cos(%s * M_PI / 180)))" % y |
97 | | - |
98 | | - def ccode_to_source(self, x, y, z): |
99 | | - return "(1000. * 1.852 * 60. * cos(%s * M_PI / 180))" % y |
100 | | - |
101 | | - |
102 | | -class GeographicSquare(UnitConverter): |
103 | | - """ Square distance converter from geometric to geographic coordinates (m2 to degree2) """ |
104 | | - source_unit = 'm2' |
105 | | - target_unit = 'degree2' |
106 | | - |
107 | | - def to_target(self, value, x, y, z): |
108 | | - return value / pow(1000. * 1.852 * 60., 2) |
109 | | - |
110 | | - def to_source(self, value, x, y, z): |
111 | | - return value * pow(1000. * 1.852 * 60., 2) |
112 | | - |
113 | | - def ccode_to_target(self, x, y, z): |
114 | | - return "pow(1.0 / (1000.0 * 1.852 * 60.0), 2)" |
115 | | - |
116 | | - def ccode_to_source(self, x, y, z): |
117 | | - return "pow((1000.0 * 1.852 * 60.0), 2)" |
118 | | - |
119 | | - |
120 | | -class GeographicPolarSquare(UnitConverter): |
121 | | - """ Square distance converter from geometric to geographic coordinates (m2 to degree2) |
122 | | - with a correction to account for narrower grid cells closer to the poles. |
123 | | - """ |
124 | | - source_unit = 'm2' |
125 | | - target_unit = 'degree2' |
126 | | - |
127 | | - def to_target(self, value, x, y, z): |
128 | | - return value / pow(1000. * 1.852 * 60. * cos(y * pi / 180), 2) |
129 | | - |
130 | | - def to_source(self, value, x, y, z): |
131 | | - return value * pow(1000. * 1.852 * 60. * cos(y * pi / 180), 2) |
132 | | - |
133 | | - def ccode_to_target(self, x, y, z): |
134 | | - return "pow(1.0 / (1000. * 1.852 * 60. * cos(%s * M_PI / 180)), 2)" % y |
135 | | - |
136 | | - def ccode_to_source(self, x, y, z): |
137 | | - return "pow((1000. * 1.852 * 60. * cos(%s * M_PI / 180)), 2)" % y |
138 | | - |
139 | | - |
140 | | -unitconverters = {'U': GeographicPolar(), 'V': Geographic(), |
141 | | - 'Kh_zonal': GeographicPolarSquare(), |
142 | | - 'Kh_meridional': GeographicSquare()} |
| 16 | +__all__ = ['Field', 'VectorField', 'SummedField', 'SummedVectorField'] |
143 | 17 |
|
144 | 18 |
|
145 | 19 | class Field(object): |
@@ -196,10 +70,10 @@ def __init__(self, name, data, lon=None, lat=None, depth=None, time=None, grid=N |
196 | 70 | self.depth = self.grid.depth |
197 | 71 | self.time = self.grid.time |
198 | 72 | fieldtype = self.name if fieldtype is None else fieldtype |
199 | | - if self.grid.mesh == 'flat' or (fieldtype not in unitconverters.keys()): |
| 73 | + if self.grid.mesh == 'flat' or (fieldtype not in unitconverters_map.keys()): |
200 | 74 | self.units = UnitConverter() |
201 | 75 | elif self.grid.mesh == 'spherical': |
202 | | - self.units = unitconverters[fieldtype] |
| 76 | + self.units = unitconverters_map[fieldtype] |
203 | 77 | else: |
204 | 78 | raise ValueError("Unsupported mesh type. Choose either: 'spherical' or 'flat'") |
205 | 79 | if type(interp_method) is dict: |
|
0 commit comments