[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[python #FTV-652122]: (No Subject)



Hi,

You can't use MetPy to import ucar.unidata.geoloc; MetPy is Python-based and 
geoloc (netCDF-java) is java-based. Those two do not play well together for the 
most part.

You *can* use MetPy together with xarray to plot/geolocate GOES ABI imagery, 
however:

import matplotlib.pyplot as plt
import metpy
from siphon.catalog import TDSCatalog

# Read data from THREDDS
cat = 
TDSCatalog('https://thredds.ucar.edu/thredds/catalog/satellite/goes/east/grb/ABI/CONUS/Channel05/current/catalog.xml')
nc = cat.datasets[0].remote_access(use_xarray=True)

# Parse the netCDF CF metadata
rad = nc.metpy.parse_cf('Rad')

# Make a plot using matplotlib and CartoPy
fig = plt.figure(figsize=(25, 15), dpi=100)
ax = plt.axes(projection=rad.metpy.cartopy_crs)
extent = (rad.x[0], rad.x[-1], rad.y[0], rad.y[-1])
ax.imshow(rad, extent=extent, cmap='Greys_r')
ax.coastlines(color='yellow', linewidth=0.5, resolution='10m')
ax.set_extent((-113, -62, 15.25, 50.75))

If you're not interested in plotting, but just getting the projection 
information out of the file, you can just use xarray and pyproj:

import matplotlib.pyplot as plt
import numpy as np
from pyproj import Proj
from siphon.catalog import TDSCatalog

# Read data from THREDDS
cat = 
TDSCatalog('https://thredds.ucar.edu/thredds/catalog/satellite/goes/east/grb/ABI/CONUS/Channel05/current/catalog.xml')
nc = cat.datasets[0].remote_access(use_xarray=True)

# Set up the PyProj projection based on metadata
h = proj_info['perspective_point_height']
geo_proj = Proj(proj='geos',
                lon_0 = proj_info['longitude_of_projection_origin'],
                lat_0 = proj_info['latitude_of_projection_origin'],
                h=h,
                a=proj_info['semi_major_axis'],
                b=proj_info['semi_minor_axis'],
                sweep=proj_info['sweep_angle_axis'],
                units='m')

# Convert 1D x, y arrays to 2D arrays. Need to scale by height to
# convert radians to meters as expected by projection
x, y = np.meshgrid(nc.Rad.x.values * h, nc.Rad.y.values * h)

# Get lat/lon values from x/y
lon, lat = geo_proj(x, y, inverse=True)

Hope this helps,

Ryan

> Dear Staff
> i wonder if I can use https://www.unidata.ucar.edu/software/metpy/
> to import ucar.unidata.geoloc
> <https://www.unidata.ucar.edu/software/netcdf-java/v4.3/v4.2/javadoc/ucar/unidata/geoloc/package-summary.html>
> libraries
> otherwise, any suggestions?
> 
> I need to do this in order to geolocate l1b abi images


Ticket Details
===================
Ticket ID: FTV-652122
Department: Support Python
Priority: Low
Status: Closed
===================
NOTE: All email exchanges with Unidata User Support are recorded in the Unidata 
inquiry tracking system and then made publicly available through the web.  If 
you do not want to have your interactions made available in this way, you must 
let us know in each email you send to us.