irfpy.util.gridsphere

A gridded sphere.

Code author: Yoshifumi Futaana

Sphere gridded into the longitude-latitude system is represented by SimpleGridSphere.

Similar to tessellatedsphere, which does not support any container functionality.

Similar to polarframe which support 3D grid and volume.

class irfpy.util.gridsphere.SimpleGridSphere(nlon=360, nlat=180)

Bases: object

A simple gridding using longitude latitude system.

A simple container in (pseudo) 2-D gridded map. Data on geometric map can be saved in it.

A grid sphere object is generated by calling SimpleGridSphere class as following.

>>> sg = SimpleGridSphere(nlon=360, nlat=180)

Longitude/Latitude pair and index

You can specify the number of grids both longitude and latitude directions.

Any data on a map can be assigned to the sphere surface. To specify the coordinates, a pair of longitude and latitude can be used. Internally, the pair is converted to a single, unique index for each grid. The index is also usedto specify the coordinate.

The conversion from longitude and latitude to the index can be done via lonlat2index().

>>> sg.lonlat2index(35.78, 28.99)
42515

This means the surface grid containing the coordinates of (lon=35.78, lat=28.99) has index of 42515.

Oppositely, the index to the center of the longitude, latitude is also available (index2lonlat()). For example, the lon and lat of surface element (center location) with id=2500 is 340.5 and -83.5, as follows.

>>> sg.index2lonlat(2500)
(340.5, -83.5)

Similar conversion from index to xyz vector is also available (index2xyz()).

Value on a map

Any data (usually scalar) can be stored by append_value_lonlat() or append_value() methods. The formar takes the pair of longitude and latitude to specify the coordinate, while the latter takes the index.

>>> sg.append_value_lonlat(25.77, -50.27, -275.25)    # Put the data -275.25 to lon=25.77 and lat=-50.27
>>> sg.append_value(2500, 275.25)

Access to the stored data

The direct access to the data by index is possible.

>>> print(sg[2500])   # Access by index
[275.25]

>>> print(sg[sg.lonlat2index(25.77, -50.27)])   # Access using lon/lat pair
[-275.25]

Indeed, any type of values can be stored; while it is recommended to use “numerical” data to store, to work with various methods.

>>> sg.append_value(3500, 'abcdefg')
>>> sg.append_value(3500, 250.5)
>>> print(sg[3500])
['abcdefg', 250.5]

The get_average() and get_counts() methods will return the 2-D array of the data. For usability, the data is already re-shaped into 2-D array. Note that the order of the index is [latitude][longitude].

>>> counts = sg.get_counts()
>>> counts.shape
(180, 360)
>>> print(counts.sum())
4

Indeed you have set a string into the data, thus, the element at index 3500 is set to np.nan when you get the average with warning message.

>>> avg = sg.get_average()
>>> print(avg.shape)
(180, 360)

Manupulations can be available get_max(), get_min() and get_median() in addition to get_average().

>>> lon, lat = sg.index2lonlat(0) 

>>> sg.append_value_lonlat(lon, lat, 30.)
>>> sg.append_value_lonlat(lon, lat, 31.)
>>> sg.append_value_lonlat(lon, lat, 32.)
>>> sg.append_value_lonlat(lon, lat, 33.)
>>> sg.append_value_lonlat(lon, lat, 34.)
>>> sg.append_value_lonlat(lon, lat, 35.)
>>> sg.append_value_lonlat(lon, lat, 43.)

Average should be 34, median to be 33, min to be 30 and max to be 43.

>>> avg = sg.get_average()
>>> print(avg[0, 0])
34.0
>>> medi = sg.get_median()
>>> print(medi[0, 0])
33.0
>>> maxx = sg.get_max()
>>> print(maxx[0, 0])
43.0
>>> minn = sg.get_min()
>>> print(minn[0, 0])
30.0

Plotting

For showing by pcolor in the matplotlib, you can just get the grid values by get_cgrid() or get_bgrid(). The first one returns the central location, but the second one returns the boundary locations.

>>> lonc, latc = sg.get_cgrid()
>>> print(lonc.shape)
(180, 360)

>>> lonb, latb = sg.get_bgrid()
>>> print(lonb.shape)
(181, 361)

To plot using pseudo-color is useful using the results of get_bgrid().

plt.pcolor(lonb, latb, avg)

Pickling

Pickling can be done.

>>> import pickle as pickle
>>> pickled = pickle.dumps(sg)
>>> sg2 = pickle.loads(pickled)
>>> sg is sg2
False
>>> sg.datalist == sg2.datalist
True
>>> sg.datalist is sg2.datalist
False

lonc

Longitude list (flattened into 1D).

latc

Latitude list (flattened into 1D).

vectorc

Central position in vector form. np.array with (3, N).

lon0

Longitude bound list

lon1

Longitude bound list

lat0

Latitude bound list

lat1

Latitude bound list

nlat

Number of latitude

nlon

Number of longitude

datalist

Datalist. Array with size of N. Each element is also array

angular_distance_matrix

A test feature of angular distance matrix. Not recommended to use. See make_angular_distance_matrix()

index2lonlat(index)

Return the (lonc, latc) for the given index.

>>> gsph = SimpleGridSphere(nlon=18, nlat=9)
>>> print(gsph.index2lonlat(120))
(250.0, 40.0)
>>> print(gsph.index2lonlat([120, 121, 122]))   
(array([250.,  270.,  290.]), array([40.,  40.,  40.]))

index2xyz(index)

lonlat2index(lon, lat)

Return the index corresponding to (lon, lat), in degrees.

>>> m = SimpleGridSphere()
>>> print(m.lonlat2index(0, -90))
0
>>> print(m.lonlat2index(1, -90))
1
>>> print(m.lonlat2index(1.5, -90))
1
>>> print(m.lonlat2index(360, -90))
0
>>> print(m.lonlat2index(360, -89))
360
>>> print(m.lonlat2index(-360, 89))
64440
>>> print(m.lonlat2index(0, 88.9))
64080
>>> print(m.lonlat2index(0, 90))
64440

lonlats2indexes(lon_list, lat_list)

get_cgrid()

Return the center grid as a pair of 2-D np.array.

Returns:

(lonarray, latarray) where lonarray (latarray) is the center longitudal (latitudal) grid in (nlat, nlon) shape.

get_bgrid(delta=0)

Return the boundary grid as a pair of 2-D np.array.

Parameters:

delta – The edges of longitude and latitude will be shrinkied with this value. If the default value 0.001 is used, the longitude range is from 0.001 to 359.999 and the latitude range is from -89.999 to 89.999. This is a bad trick, but s very useful for plotting.

Returns:

(lonarray, latarray) where lonarray (latarray) is the boundary longitudal (latitudal) grid in (nlat+1, nlon+1) shape.

get_average()

Get average of each grid.

Returns:

np.ma.masked_array with shape of (nlat, nlon)

Return type:

np.ma.masked_array

It works correctly only when all the data is float (or equivalent). Otherwise, masked.

>>> gs = SimpleGridSphere()
>>> gs.append_value_lonlat(150, 65, 25.3)
>>> gs.append_value_lonlat(358, 27.5, 883)
>>> gs.append_value_lonlat(25, 90, -35)

>>> lons, lats = gs.get_bgrid()

>>> ave = gs.get_average()

You can use pcolor with:

plt.pcolor(lons, lats, ave)

get_average2()

Return the averaged value of the data in each grid.

It is equivalent to get_average(), but more primitive implementaiton. Recommended to use get_average() for usual purpose.

get_median()

get_max()

get_min()

get_counts()

Return the number of data stored in the grid

Similar to get_average(), while usual np.array is returned.

append_value(index, value)

Append the given value to the given index.

append_values(indexlist, valuelist)

append_value_lonlat(lon, lat, value)

Append the given value to the given (lon, lat).

Parameters:

• lon – Longitude (degrees)

• lat – Latitude (degrees)

• value – Value

append_values_lonlat(lonlist, latlist, valuelist)

clear()

Clear the data. All the list in the data will be replaced by [].

make_angular_distance_matrix(sqrt=True, clear_cache=False)

Create a distance matrix.

The distance matrix is the matrix that saves the distance between the elements. It can be used, for example, to calculate the spacecraft visibility from the surface of body.

Warning

The method will consume time and memory quite a lot.

This functionality is memory hungry. Indeed, it consumes (nlon * nlat) ** 2. This means that, as default gridding (360, 180) will need 4.2e9 data, which may consume more than 10 GB memory!! However, for smaller gridding, e.g. (36, 18) every 10 degrees, may reduces 4.2e5 data, which may consume only a few MB.

This functionality also time consuming. On iMac in my office it takes ~1 sec for nlon*nlat = 144. The algorithm is N(O^2), thus, (nlon=36, nlat=18) will take ~25 sec.

Generated matrix will have (N, N) elements, where N is nlon * nlat. matrix[i, j] have the minimum distance between the grid indexed i and j for i < j. matrix[i, j] have the maximum distance between the grid indexed i and j for i > j. matrix[i, j] have zero for i == j.

>>> gs = SimpleGridSphere(nlon=6, nlat=3)
>>> gs.make_angular_distance_matrix()