open_nc.txt

import xarray as xr
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import numpy as np

# Load the NetCDF file
file_path = r'file_directory'
data = xr.open_dataset(file_path)

# Choose the variable to visualize
variable_name = 'VODCA_L'
data_variable = data[variable_name]

# Check if the 'time' dimension exists and select the first time slice if it does
if 'time' in data_variable.dims:
    data_variable = data_variable.isel(time=0)

# Create a plot with a world map background
plt.figure(figsize=(12, 6))
ax = plt.axes(projection=ccrs.PlateCarree())

# Plot the data
data_variable.plot(ax=ax, transform=ccrs.PlateCarree(), cmap='viridis', cbar_kwargs={'label': 'Vegetation Optical Depth'})

# Add coastlines and other features using built-in Cartopy features
ax.coastlines()
ax.add_feature(cfeature.BORDERS, linestyle=':')
ax.add_feature(cfeature.LAND, facecolor='lightgray')
ax.add_feature(cfeature.OCEAN, facecolor='lightblue')
ax.add_feature(cfeature.LAKES, edgecolor='black', facecolor='lightblue')
ax.add_feature(cfeature.RIVERS)

# Convert the numpy.datetime64 object to a datetime object if time exists
if 'time' in data:
    date_str = np.datetime_as_string(data['time'].values[0], unit='D')
    plt.title(f"{variable_name} on {date_str}")
else:
    plt.title(variable_name)

plt.xlabel('Longitude')
plt.ylabel('Latitude')
plt.show()