1

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)
Marshall Sutton
  • 21
  • 2
  • Have you tried the transform [that worked for this question](https://stackoverflow.com/questions/60095078/cartopy-coastlines-and-contourf-interfering?rq=1)? I would also try changing the z-order of whichever you want in back, or whichever takes z-order, depending. I would try it out myself on this very tidy example but I'm having a hassle installing cartopy. – cphlewis Jan 29 '22 at 21:07
  • @cphlewis You need to install GEOS before you install Cartopy. https://libgeos.org/ Adding the transform from the question helped a little. I can see a smattering of data plotted over the coastlines. – Marshall Sutton Jan 30 '22 at 17:59
  • @Marshall Sutton, you still working on this? If so, please clarify how you get "_can either plot the coastlines, or the data, but not both_". I can get the coastlines only, but data doesn't display with `ax.pcolormesh()`. How did you plot the data? I am able to plot the image data with `plt.imshow(cloud_data)`. The display is better if you mask the cloud data using `_FillValue`. – kcw78 Mar 12 '22 at 17:25
  • @kcw78 I used the transform that worked in the comment above. I also had to zero all the non-earth pixels. They were set to a value of -999. Setting them to a value of zero seemed to fix the problem. Here is my latest solution. – Marshall Sutton Mar 16 '22 at 17:15

1 Answers1

1

@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)
Marshall Sutton
  • 21
  • 2