3

I have a large netCDF file (444GB) that I would like to transpose. The file is flipped with latitude on the x-axis and lon on the y-axis. See pictures of one band below.

not correct

Not correct

Correct

Correct: I transposed this using raster package in R but my 444GB file is too big and the process has been running for 30 hours now non-stop. Thus, it is too slow to accomplish this task. How can I achieve the same task using GDAL, NCL, NCO, or CDO?

An example file can be gotten from here link to files

2

You mention using R - would a Python solution be suitable? I'm not sure how this will compare in terms of processing time.

The netCDF4 Python library allows you to store your variable in a numpy array. You can then use numpy to manipulate the array, in your case you can use the numpy.transpose() method.

The following shows how to access a 2D array from a netCDF file and save it as a GeoTIFF file:

import gdal,gdalconst,osr
import numpy
from netCDF4 import Dataset

# Open the NetCDF file for reading.
fh = Dataset(r"INPUT_PATH.nc","r")

# Read your parameter into a numpy array.
# Depending on your dataset you may have to specify additional index values
# e.g. Time.
values = numpy.array(fh.variables["YOUR_PARAMETER"])

# Now you can manipulate your array.
"""
Transpose method goes here
"""

# Get array dimensions so we can create an equal size tiff file.
rows,cols = values.shape

# Create output tiff file.
# driver.Create() parameters are: output path, number of columns, number of rows,
# number of bands, data type
driver = gdal.GetDriverByName("GTiff")
out_data = r"OUTPUT_PATH.tif"
out_tif = driver.Create(out_data, cols, rows, 1, gdalconst.GDT_Float32)

# Create Spatial Reference object and set GeoTIFF projection.
# This information may be found in either the data documentation or the netCDF file.
prj = osr.SpatialReference()
prj.ImportFromEPSG(4326) # WGS84
out_tif.SetProjection(prj.ExportToWkt())

# Set GeoTransformation.
# This information may be found in either the data documentation, the netCDF file, or
# can be derived. For example, if you know the longitude range and number of columns
# you can calculate the x step as float(lon_range)/float(num_cols).
geotrans = [top_left_x,x_step,rotation,top_left_y,rotation,y_step]
out_tif.SetGeoTransform(geotrans)

# Finally we can write the array to the raster band.
out_band = out_tif.GetRasterBand(1)
out_band.WriteArray(values)

# Clear the memory and close the output file.
out_tif.FlushCache()
out_tif = None

Your main analysis will go in the "Transpose method goes here" section. When transposing the array remember that you will need to account for these changes in the output GeoTIFF creation and geotransformation. For example the number of rows and columns will be reversed, as will the starting coordinates and x and y steps.

  • thanks for providing a more robust solution within python. Unfortunately I have no working knowledge of python. I use mostly GDAL, NCL, NCO, CDO, MATLAB and R. – code123 Mar 17 '16 at 16:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.