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

[python #FDY-952600]: Building Matplotlib Images using AWIPS Python Framework

Hi,

Apologies on the delayed response, I've been buried. So the problem that I see 
is that if you're going to plot using imshow, you have to have the native 
projection coordinates for the grid of data. That's because all you give imshow 
for positioning is the bounding box of the image, so it assumes each pixel is 
uniform in size, evenly dividing the rectangle you pass it with extent. 
Fortunately, using the Geostationary CRS you set up with CartoPy, you can 
create those x,y values from the lons/lats you get from EDEX. Here's an example 
that seems to work for me:

from awips.dataaccess import DataAccessLayer
import cartopy.feature as cfeature
import cartopy.crs as ccrs
import matplotlib.pyplot as plt
import numpy as np

edexServer='edex-cloud.unidata.ucar.edu'
DataAccessLayer.changeEDEXHost(edexServer)
request = DataAccessLayer.newDataRequest()
request.setDatatype("satellite")
sector="ECONUS"
request.setLocationNames(sector)
request.setParameters("CH-13-10.35um")
cycles = DataAccessLayer.getAvailableTimes(request, True)
times = DataAccessLayer.getAvailableTimes(request)
fcstRun = DataAccessLayer.getForecastRun(cycles[-1], times)
response = DataAccessLayer.getGridData(request, [fcstRun[0]])
grid = response[0]
imdata = grid.getRawData()

# Set up the projection using known parameters
globe = ccrs.Globe(semimajor_axis=6378137.0, semiminor_axis=6356752.31414)
crs = ccrs.Geostationary(central_longitude=-75,
                         satellite_height=35786023, globe=globe, sweep_axis='x')

# Convert grid lons/lats to native projected coordinates
lons, lats = grid.getLatLonCoords()
x, y, _ = crs.transform_points(ccrs.PlateCarree(), lons, lats).T

bbox = [-130.0, -65.0, 20.0, 55.0]
fig = plt.figure(figsize=(6,6),dpi=200)
ax = plt.axes(projection=ccrs.PlateCarree())
ax.set_extent(bbox, ccrs.PlateCarree())
ax.add_feature(cfeature.COASTLINE)
ax.add_feature(cfeature.STATES)

image = ax.imshow(imdata, origin='upper', extent=(np.nanmin(x), np.nanmax(x), 
np.nanmin(y), np.nanmax(y)), cmap='Greys', vmin=170, vmax=375, transform=crs)

Note I also added the "sweep_axis='x'" parameter to the construction of 
Geostationary as well as tweaked the value of `satellite_height`. Everything 
else, as far as masked values and whatnot does not appear to be a problem (note 
the use of nanmin/nanmax to avoid manual masking).

Hope this helps,

Ryan

> Hi Shay and Ryan,
> 
> I was wondering if one of you would be able to assist me with a python
> script to build "imshow" images from SBN GOES-16 data sets?
> 
> I am using your awips example provided in the tutorial for satellite
> imagery,
> <https://python-awips.readthedocs.io/en/latest/examples/generated/Satellite_Imagery.html>
> but am running into an issue with plotting non-finite data and plotting the
> data into the correct projection. Do you know what projection the data is
> disseminated in on SBN or how I could resolve this?
> 
> My goal is to capture an image using imshow... mainly because pcolormesh
> takes quite a long time to map images when dealing with large shapes and
> the nans. Attached shows the correct image made with metpy and the
> incorrect image I am getting with the awips python framework.
> 
> Any help would be greatly appreciated.

Ticket Details
===================
Ticket ID: FDY-952600
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.