Question

I am facing difficulties with reprojecting a NetCDF file containing data on a rotated latitude-longitude grid. Despite several attempts, I have been unable to align the data correctly with geographic locations (such as Italy) when visualizing in QGIS. I have tried various GDAL commands and manual approaches without success.

NC FILE LINK

Data summary:

Dimensions:                         (time: 24, rlat: 652, rlon: 570, bnds: 2)
Coordinates:
  * time                            (time) datetime64[ns] 192B 2006-01-01 ......
  * rlon                            (rlon) float64 5kB -5.913 -5.893 ... 5.467
  * rlat                            (rlat) float64 5kB -6.408 -6.388 ... 6.612
    lat                             (rlat, rlon) float32 1MB ...
    lon                             (rlat, rlon) float32 1MB ...
Dimensions without coordinates: bnds
Data variables:
    TOT_PREC                        (time, rlat, rlon) float32 36MB ...
    crs_rotated_latitude_longitude  int64 8B ...
    time_bnds                       (time, bnds) datetime64[ns] 384B ...
Attributes: (12/13)
    CDI:             Climate Data Interface version 1.9.8 (https://mpimet.mpg...
    Conventions:     CF-1.4
    history:         Generated by CMCC DDS version 0.9.0 2024-07-29 08:33:07....
    source:          COSMO
    institution:     CMCC (Euro-Mediterranean Center on Climate Change), REMH...
    title:           cclm-sp_2.4_terra_urb_2.3.1 simulation (0.02 Deg) forced...
    ...              ...
    experiment_id:   Historical with Urban Parametrization, DT=20sec
    conventionsURL:  http://www.cfconventions.org/
    contact:         Mario Raffa ([email protected])
    references:      http//www.clm-community.eu, http://www.cmcc.it
    creation_date:   2021-12-31 12:05:48
    CDO:             Climate Data Operators version 1.9.8 (https://mpimet.mpg...

Commands attempted:

Using gdal_translate with SRS:

gdal_translate -a_srs "EPSG:4326" precipitation.nc precipitation_with_srs.nc

Error: NetCDF: Index exceeds dimension bound.

I have also tried to handle the transformation using Python with xarray, numpy, and pyproj. Here’s the code I used to transform and save the NetCDF file:

import xarray as xr
import numpy as np
from pyproj import Transformer
import os

# Path to the original NetCDF file
input_file_path = r'C:UsersHUNDESADesktophun - Copy
precipitation_amount_era5-downscaled-over-italy_hourly_418082.nc'

# Path for the new NetCDF file
output_file_path = r'C:UsersHUNDESADesktophun - Copytransformed_file.nc'

# Open the original NetCDF file using xarray
ds = xr.open_dataset(input_file_path)

# Rotation parameters
north_pole_lat = 47
north_pole_lon = -168

# Define the coordinate system transformation
transformer = Transformer.from_proj(
    proj_from=f"+proj=aeqd +lat_0={north_pole_lat} +lon_0={north_pole_lon}",
    proj_to="epsg:4326"  # WGS84
)

# Create a copy of the dataset
new_ds = ds.copy()

# Transform the coordinates if necessary


# Save the new dataset to a new NetCDF file
encoding = {
    'TOT_PREC': {'dtype': 'float32'},
    'latitude': {'dtype': 'float32'},
    'longitude': {'dtype': 'float32'},
}

# Convert DataArray to DataArray with proper dtype if needed
for var_name in encoding.keys():
    if var_name in new_ds.variables:
        new_ds[var_name] = new_ds[var_name].astype(encoding[var_name]['dtype'])

# Save the dataset
new_ds.to_netcdf(output_file_path, encoding=encoding)

Error: I’m unable to properly reproject the data to EPSG:4326, resulting in distorted visualization.

Trouble Reprojecting Rotated Latitude-Longitude Grid to EPSG:4326 in GDAL

LEAVE A COMMENT Hủy