|
| 1 | +""" |
| 2 | +project - Project data onto lines or great circles, or generate tracks. |
| 3 | +""" |
| 4 | +import pandas as pd |
| 5 | +from pygmt.clib import Session |
| 6 | +from pygmt.exceptions import GMTInvalidInput |
| 7 | +from pygmt.helpers import ( |
| 8 | + GMTTempFile, |
| 9 | + build_arg_string, |
| 10 | + fmt_docstring, |
| 11 | + kwargs_to_strings, |
| 12 | + use_alias, |
| 13 | +) |
| 14 | + |
| 15 | + |
| 16 | +@fmt_docstring |
| 17 | +@use_alias( |
| 18 | + A="azimuth", |
| 19 | + C="center", |
| 20 | + E="endpoint", |
| 21 | + F="convention", |
| 22 | + G="generate", |
| 23 | + L="length", |
| 24 | + N="flat_earth", |
| 25 | + Q="unit", |
| 26 | + S="sort", |
| 27 | + T="pole", |
| 28 | + V="verbose", |
| 29 | + W="width", |
| 30 | + Z="ellipse", |
| 31 | + f="coltypes", |
| 32 | +) |
| 33 | +@kwargs_to_strings(E="sequence", L="sequence", T="sequence", W="sequence", C="sequence") |
| 34 | +def project(data=None, x=None, y=None, z=None, outfile=None, **kwargs): |
| 35 | + r""" |
| 36 | + Project data onto lines or great circles, or generate tracks. |
| 37 | +
|
| 38 | + Project reads arbitrary :math:`(x, y [, z])` data and returns any |
| 39 | + combination of :math:`(x, y, z, p, q, r, s)`, where :math:`(p, q)` are the |
| 40 | + coordinates in the projection, :math:`(r, s)` is the position in the |
| 41 | + :math:`(x, y)` coordinate system of the point on the profile (:math:`q = 0` |
| 42 | + path) closest to :math:`(x, y)`, and :math:`z` is all remaining columns in |
| 43 | + the input (beyond the required :math:`x` and :math:`y` columns). |
| 44 | +
|
| 45 | + Alternatively, ``project`` may be used to generate |
| 46 | + :math:`(r, s, p)` triples at equal increments along a profile using the |
| 47 | + ``generate`` parameter. In this case, the value of ``data`` is ignored |
| 48 | + (you can use, e.g., ``data=None``). |
| 49 | +
|
| 50 | + Projections are defined in any (but only) one of three ways: |
| 51 | +
|
| 52 | + 1. By a ``center`` and an ``azimuth`` in degrees clockwise from North. |
| 53 | + 2. By a ``center`` and ``endpoint`` of the projection path. |
| 54 | + 3. By a ``center`` and a ``pole`` position. |
| 55 | +
|
| 56 | + To spherically project data along a great circle path, an oblique |
| 57 | + coordinate system is created which has its equator along that path, and the |
| 58 | + zero meridian through the Center. Then the oblique longitude (:math:`p`) |
| 59 | + corresponds to the distance from the Center along the great circle, and the |
| 60 | + oblique latitude (:math:`q`) corresponds to the distance perpendicular to |
| 61 | + the great circle path. When moving in the increasing (:math:`p`) direction, |
| 62 | + (toward B or in the azimuth direction), the positive (:math:`q`) direction |
| 63 | + is to your left. If a Pole has been specified, then the positive |
| 64 | + (:math:`q`) direction is toward the pole. |
| 65 | +
|
| 66 | + To specify an oblique projection, use the ``pole`` option to set |
| 67 | + the pole. Then the equator of the projection is already determined and the |
| 68 | + ``center`` option is used to locate the :math:`p = 0` meridian. The center |
| 69 | + *cx/cy* will be taken as a point through which the :math:`p = 0` meridian |
| 70 | + passes. If you do not care to choose a particular point, use the South pole |
| 71 | + (*cx* = 0, *cy* = -90). |
| 72 | +
|
| 73 | + Data can be selectively windowed by using the ``length`` and ``width`` |
| 74 | + options. If ``width`` is used, the projection width is set to use only |
| 75 | + data with :math:`w_{{min}} < q < w_{{max}}`. If ``length`` is set, then |
| 76 | + the length is set to use only those data with |
| 77 | + :math:`l_{{min}} < p < l_{{max}}`. If the ``endpoint`` option |
| 78 | + has been used to define the projection, then ``length="w"`` may be used to |
| 79 | + window the length of the projection to exactly the span from O to B. |
| 80 | +
|
| 81 | + Flat Earth (Cartesian) coordinate transformations can also be made. Set |
| 82 | + ``flat_earth=True`` and remember that azimuth is clockwise from North (the |
| 83 | + y axis), NOT the usual cartesian theta, which is counterclockwise from the |
| 84 | + x axis. azimuth = 90 - theta. |
| 85 | +
|
| 86 | + No assumptions are made regarding the units for |
| 87 | + :math:`x, y, r, s, p, q, dist, l_{{min}}, l_{{max}}, w_{{min}}, w_{{max}}`. |
| 88 | + If -Q is selected, map units are assumed and :math:`x, y, r, s` must be in |
| 89 | + degrees and :math:`p, q, dist, l_{{min}}, l_{{max}}, w_{{min}}, w_{{max}}` |
| 90 | + will be in km. |
| 91 | +
|
| 92 | + Calculations of specific great-circle and geodesic distances or for |
| 93 | + back-azimuths or azimuths are better done using :gmt-docs:`mapproject` as |
| 94 | + project is strictly spherical. |
| 95 | +
|
| 96 | + {aliases} |
| 97 | +
|
| 98 | + Parameters |
| 99 | + ---------- |
| 100 | + data : str or {table-like} |
| 101 | + Pass in (x, y, z) or (longitude, latitude, elevation) values by |
| 102 | + providing a file name to an ASCII data table, a 2D |
| 103 | + {table-classes}. |
| 104 | +
|
| 105 | + center : str or list |
| 106 | + *cx*/*cy*. |
| 107 | + Set the origin of the projection, in Definition 1 or 2. If |
| 108 | + Definition 3 is used, then *cx/cy* are the coordinates of a |
| 109 | + point through which the oblique zero meridian (:math:`p = 0`) should |
| 110 | + pass. The *cx/cy* is not required to be 90 degrees from the pole. |
| 111 | +
|
| 112 | + azimuth : float or str |
| 113 | + Define the azimuth of the projection (Definition 1). |
| 114 | +
|
| 115 | + endpoint : str or list |
| 116 | + *bx*/*by*. |
| 117 | + Define the end point of the projection path (Definition 2). |
| 118 | +
|
| 119 | + convention : str |
| 120 | + Specify the desired output using any combination of **xyzpqrs**, in |
| 121 | + any order [Default is **xypqrsz**]. Do not space between the letters. |
| 122 | + Use lower case. The output will be columns of values corresponding to |
| 123 | + your ``convention``. The **z** flag is special and refers to all |
| 124 | + numerical columns beyond the leading **x** and **y** in your input |
| 125 | + record. The **z** flag also includes any trailing text (which is |
| 126 | + placed at the end of the record regardless of the order of **z** in |
| 127 | + ``convention``). **Note**: If ``generate`` is True, then the output |
| 128 | + order is hardwired to be **rsp** and ``convention`` is not allowed. |
| 129 | +
|
| 130 | + generate : str |
| 131 | + *dist* [/*colat*][**+c**\|\ **h**]. |
| 132 | + Create :math:`(r, s, p)` output data every *dist* units of :math:`p` |
| 133 | + (See `unit` option). Alternatively, append */colat* for a small |
| 134 | + circle instead [Default is a colatitude of 90, i.e., a great circle]. |
| 135 | + If setting a pole with ``pole`` and you want the small circle to go |
| 136 | + through *cx*/*cy*, append **+c** to compute the required colatitude. |
| 137 | + Use ``center`` and ``endpoint`` to generate a circle that goes |
| 138 | + through the center and end point. Note, in this case the center and |
| 139 | + end point cannot be farther apart than :math:`2|\mbox{{colat}}|`. |
| 140 | + Finally, if you append **+h** then we will report the position of |
| 141 | + the pole as part of the segment header [Default is no header]. |
| 142 | + Note: No input is read and the value of ``data``, ``x``, ``y``, |
| 143 | + and ``z`` is ignored if ``generate`` is used. |
| 144 | +
|
| 145 | + length : str or list |
| 146 | + [**w**\|\ *l_min*/*l_max*]. |
| 147 | + Project only those data whose *p* coordinate is |
| 148 | + within :math:`l_{{min}} < p < l_{{max}}`. If ``endpoint`` has been set, |
| 149 | + then you may alternatively use **w** to stay within the distance from |
| 150 | + ``center`` to ``endpoint``. |
| 151 | +
|
| 152 | + flat_earth : bool |
| 153 | + Make a Cartesian coordinate transformation in the plane. |
| 154 | + [Default is ``False``; plane created with spherical trigonometry.] |
| 155 | +
|
| 156 | + unit : bool |
| 157 | + Set units for :math:`x, y, r, s` degrees and |
| 158 | + :math:`p, q, dist, l_{{min}}, l_{{max}}, w_{{min}}, {{w_max}}` to km. |
| 159 | + [Default is ``False``; all arguments use the same units] |
| 160 | +
|
| 161 | + sort : bool |
| 162 | + Sort the output into increasing :math:`p` order. Useful when projecting |
| 163 | + random data into a sequential profile. |
| 164 | +
|
| 165 | + pole : str or list |
| 166 | + *px*/*py*. |
| 167 | + Set the position of the rotation pole of the projection. |
| 168 | + (Definition 3). |
| 169 | +
|
| 170 | + {V} |
| 171 | +
|
| 172 | + width : str or list |
| 173 | + *w_min*/*w_max*. |
| 174 | + Project only those data whose :math:`q` coordinate is |
| 175 | + within :math:`w_{{min}} < q < w_{{max}}`. |
| 176 | +
|
| 177 | + ellipse : str |
| 178 | + *major*/*minor*/*azimuth* [**+e**\|\ **n**]. |
| 179 | + Used in conjunction with ``center`` (sets its center) and ``generate`` |
| 180 | + (sets the distance increment) to create the coordinates of an ellipse |
| 181 | + with *major* and *minor* axes given in km (unless ``flat_earth`` is |
| 182 | + given for a Cartesian ellipse) and the *azimuth* of the major axis in |
| 183 | + degrees. Append **+e** to adjust the increment set via ``generate`` so |
| 184 | + that the the ellipse has equal distance increments [Default uses the |
| 185 | + given increment and closes the ellipse]. Instead, append **+n** to set |
| 186 | + a specific number of unique equidistant data via ``generate``. For |
| 187 | + degenerate ellipses you can just supply a single *diameter* instead. A |
| 188 | + geographic diameter may be specified in any desired unit other than km |
| 189 | + by appending the unit (e.g., 3d for degrees) [Default is km]; |
| 190 | + the increment is assumed to be in the same unit. **Note**: |
| 191 | + For the Cartesian ellipse (which requires ``flat_earth``), the |
| 192 | + *direction* is counter-clockwise from the horizontal instead of an |
| 193 | + *azimuth*. |
| 194 | +
|
| 195 | + outfile : str |
| 196 | + The file name for the output ASCII file. |
| 197 | +
|
| 198 | + {f} |
| 199 | +
|
| 200 | + Returns |
| 201 | + ------- |
| 202 | + track: pandas.DataFrame or None |
| 203 | + Return type depends on whether the ``outfile`` parameter is set: |
| 204 | +
|
| 205 | + - :class:`pandas.DataFrame` table with (x, y, ..., newcolname) if |
| 206 | + ``outfile`` is not set |
| 207 | + - None if ``outfile`` is set (output will be stored in file set |
| 208 | + by ``outfile``) |
| 209 | + """ |
| 210 | + |
| 211 | + if "C" not in kwargs: |
| 212 | + raise GMTInvalidInput("The `center` parameter must be specified.") |
| 213 | + if "G" not in kwargs and data is None: |
| 214 | + raise GMTInvalidInput( |
| 215 | + "The `data` parameter must be specified unless `generate` is used." |
| 216 | + ) |
| 217 | + if "G" in kwargs and "F" in kwargs: |
| 218 | + raise GMTInvalidInput( |
| 219 | + "The `convention` parameter is not allowed with `generate`." |
| 220 | + ) |
| 221 | + |
| 222 | + with GMTTempFile(suffix=".csv") as tmpfile: |
| 223 | + if outfile is None: # Output to tmpfile if outfile is not set |
| 224 | + outfile = tmpfile.name |
| 225 | + with Session() as lib: |
| 226 | + if "G" not in kwargs: |
| 227 | + # Choose how data will be passed into the module |
| 228 | + table_context = lib.virtualfile_from_data( |
| 229 | + check_kind="vector", data=data, x=x, y=y, z=z, required_z=False |
| 230 | + ) |
| 231 | + |
| 232 | + # Run project on the temporary (csv) data table |
| 233 | + with table_context as infile: |
| 234 | + arg_str = " ".join( |
| 235 | + [infile, build_arg_string(kwargs), "->" + outfile] |
| 236 | + ) |
| 237 | + else: |
| 238 | + arg_str = " ".join([build_arg_string(kwargs), "->" + outfile]) |
| 239 | + lib.call_module(module="project", args=arg_str) |
| 240 | + |
| 241 | + # if user did not set outfile, return pd.DataFrame |
| 242 | + if outfile == tmpfile.name: |
| 243 | + if "G" in kwargs: |
| 244 | + column_names = list("rsp") |
| 245 | + result = pd.read_csv(tmpfile.name, sep="\t", names=column_names) |
| 246 | + else: |
| 247 | + result = pd.read_csv(tmpfile.name, sep="\t", header=None, comment=">") |
| 248 | + # return None if outfile set, output in outfile |
| 249 | + elif outfile != tmpfile.name: |
| 250 | + result = None |
| 251 | + |
| 252 | + return result |
0 commit comments