Marshall Sutton
Marshall Sutton

Reputation: 21

Plotting Satellite data and coastlines using Matplotlib and Cartopy

I am trying to plot data from the DSCOVR Satellite onto an orthographic projection and add coastlines to the image. I am using Matplotlib and and cartopy. I can either get the coastlines to show up, or the data, but not both. I suspect that part of the problem is that the Latitude and Longitude data include dark areas (the space that surrounds the earth) and are given the value -999.0.

Data can be found at: https://asdc.larc.nasa.gov/data/DSCOVR/EPIC/L2_CLOUD_03/2021/12/ Any file from that directory will work.

Here is my code:

import cartopy.crs as ccrs
import h5py
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap
import numpy as np


filename = 'DSCOVR_EPIC_L2_CLOUD_03_20211203163932_03.nc4'

def make_cloud_image(cloud_name):
    cloud = h5py.File(filename, 'r')
    cloud_data = np.array(cloud['geophysical_data/Cloud_Mask'])
    lat_data = np.array(cloud['geolocation_data/latitude'])
    lon_data = np.array(cloud['geolocation_data/latitude'])
    lat = cloud.attrs.__getitem__('centroid_mean_latitude')[0]
    lon = cloud.attrs.__getitem__('centroid_mean_longitude')[0]
    projection = ccrs.Orthographic(central_longitude=lon, central_latitude=lat, globe=None)    
    ax = plt.subplot(projection=projection)
    cmap = ListedColormap(["Black","dodgerblue", "lightblue", "lightyellow", "white"])
    ax.pcolormesh(lat_data,lon_data, cloud_data, cmap=cmap, transform=projection)
    ax.annotate(cloud_name[-30:-4],xy=(.5,.90),xycoords='axes fraction',
                color='white',fontsize=6, horizontalalignment='center')
    ax.set_global()
    plt.axis('off')
    ax.coastlines()

    plt.show
    

make_cloud_image(filename)

Upvotes: 1

Views: 1392

Answers (1)

Marshall Sutton
Marshall Sutton

Reputation: 21

@kcw78, I had to zero out the non-earth pixels, as they were set to a value of '-999' I also used the transform suggested in the earlier comment.

#import cartopy.crs as ccrs

import h5py
import os, sys
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap
from matplotlib.patches import Patch
import numpy as np
import cartopy.crs as ccrs

filename = 'DSCOVR_EPIC_L2_CLOUD_03_20211203163932_03.nc4'
fontsiz = 6



def make_cloud_image(cloud_name):
    cloud = h5py.File(filename, 'r')
    data = np.array(cloud['geophysical_data/Cloud_Mask'])
    data_where = np.where(data==0,np.nan,data)
    lat_data = np.array(cloud['geolocation_data/latitude'])
    lat_where = np.where(lat_data==-999.0,0,lat_data)
    lon_data = np.array(cloud['geolocation_data/longitude'])
    lon_where = np.where(lon_data==-999.0,0,lon_data)
    lat = cloud.attrs.__getitem__('centroid_mean_latitude')[0]
    lon = cloud.attrs.__getitem__('centroid_mean_longitude')[0]
    lon_min = cloud.attrs.__getitem__('minimum_longitude')[0]
    lat_min = cloud.attrs.__getitem__('minimum_latitude')[0]    
    lon_max = cloud.attrs.__getitem__('maximum_longitude')[0]
    lat_max = cloud.attrs.__getitem__('maximum_latitude')[0]
    projection=ccrs.Orthographic(central_longitude=lon, central_latitude=lat,globe=None)
    
    ax = plt.subplot(projection=projection)
        
    cmap = ListedColormap(["dodgerblue", "lightblue", "lightyellow", "white"])
    plt.pcolormesh(lon_where,lat_where,data_where, cmap=cmap, transform=ccrs.PlateCarree(), zorder=2)
                     
    ax.coastlines(zorder=3)
    
    plt.savefig(f'{filename[:-4]}.png')
    
make_cloud_image(filename)

Upvotes: 1

Related Questions