1

I'm working with UKCP18 climate data (2.2km), which uses an unfamiliar georeferencing (see link https://www.metoffice.gov.uk/pub/data/weather/uk/ukcp18/science-reports/UKCP-Convection-permitting-model-projections-report.pdf at Section 2.3 for desc.), and am trying to convert the file into a TIFF to use, but with traditional lat, lon (aka, EPSG:4326).

It uses two spatial dimensions "grid_latitude" and "grid_longitude" (in this unfamiliar georeferenced grid) and I'm looking to convert the file into the non-dimensional coordinates "latitude" and "longitude" present within the file, both of which are functions of (grid_latitude, grid_longitude). I'm new to working with files of this type, but want to avoid unnecessary reprojections later in gdalwarp if the file itself actually contains lat, lon info like this one already does; as 2D coordinates.

How do I find a way to set the spatial dimensions as latitude and longitude when converting to a raster? It seems lots of people online have had this issue, but I have not seen a successful fix to the issue.

Here is the current state of the xr.DataArray after working on the file:

<xarray.DataArray 'tasmin' (grid_latitude: 606, grid_longitude: 484)>
dask.array<mean_agg-aggregate, shape=(606, 484), dtype=float32, chunksize=(606, 484), chunktype=numpy.ndarray>
Coordinates:
  * grid_latitude   (grid_latitude) float64 -4.683 -4.647 -4.611 ... 8.027 8.063
  * grid_longitude  (grid_longitude) float64 353.9 354.0 354.0 ... 364.3 364.3
    latitude        (grid_latitude, grid_longitude) float64 dask.array<chunksize=(606, 484), meta=np.ndarray>
    longitude       (grid_latitude, grid_longitude) float64 dask.array<chunksize=(606, 484), meta=np.ndarray>

If I can get the file in the format:

<xarray.DataArray 'tasmin' (latitude: 606, longitude: 484)>

... that would be perfect.

The file itself recognises how this would be georeferenced with lat, lon as I can plot it with:

da.plot(x='longitude', y='latitude')

which gives: image of gridded minimum temperatures over the United Kingdom on an irregular mesh

Thanks in advance!

Michael Delgado
  • 13,789
  • 3
  • 29
  • 54
Kieran Sam Bajpai
  • 21
  • 2
  • this looks like a job for [`xESMF`](https://xesmf.readthedocs.io/en/latest/index.html)! Welcome to the wild world of projections and non-equirectangular grids :) – Michael Delgado Jun 02 '22 at 19:19
  • Cheers @MichaelDelgado - it looks useful! I am just confused that there's not a 1/2 step process to just say "can I save this file but with lat,lon instead? lat,lon are literally coordinates in the file? They're just 2D – Kieran Sam Bajpai Jun 02 '22 at 20:17
  • unfortunately not. your data is at regular intervals in (x, y), and after a highly nonlinear transform, these x, y points correspond to lat/lon values which are not at all regularly spaced. if you want to force the data to be on a regular lat/lon grid (which is a totally normal thing to do) you need to decide how you want to interpolate to your desired regular grid. there will have to be some resampling/interpolating/data destruction in this process, as you're imposing a new grid specification and will have to use some method to get values at each of the new points. – Michael Delgado Jun 02 '22 at 20:33
  • see this figure of [an irregular mesh](//xesmf.readthedocs.io/en/latest/notebooks/Curvilinear_grid.html) in the xesmf docs. each of the grid cells on this mesh is defined by an x, y point, but when projected onto an equirectangular lat/lon plot like EPSG:4326 you can see how the x/y to lat/lon conversion is not at all just a matter of swapping out the dims. and in that example, you can see how the choice of an output grid and sampling algorithm really does matter - the example outputs are very blocky and a lot of data has been lost. but this may be the right choice depending on your use case! – Michael Delgado Jun 02 '22 at 20:37
  • the reason you can't just save the file with dims (lat, lon) instead of (x, y), or why you can't have 2D dims in general, is because the dims are the dimensions of the actual array in memory/on disk. You can't have a numpy array of shape (lat(x, y), lon(x, y)). It's just shape (len(y), len(x)), unless you regrid the data so it's shape (len(lat), len(lon)). – Michael Delgado Jun 02 '22 at 20:42
  • 1
    Thanks for all the information and the modifications @MichaelDelgado! Very much appreciated and I'm now experimenting between xESMF and gdalwarp for results and quality! – Kieran Sam Bajpai Jun 03 '22 at 14:32
  • It all worked brilliantly when I used xESMF @MichaelDelgado - there were a few hurdles afterwards to get over, like recognising how to use the xe.util.... function and creating a new xr.Dataset, but I got there. You've saved me countless hours, so thanks again. – Kieran Sam Bajpai Jun 03 '22 at 18:49
  • great! if you're up for it - you could post your solution here - I'm sure it would be a valuable guide for others :) – Michael Delgado Jun 03 '22 at 20:53

0 Answers0