I am producing a raster file in python. However, the rasters are incorrect: they appear to be rotated 90° clockwise. I have provided sample code and results to illustrate the problem.
This is basically how I wanted my raster to look like (using his2d).
But the raster, when opened in Qgis, looks like this:
As you can see, it seems to be rotated 90 degrees clockwise. Below is an example of the code I am using. Note I have just created random data because the actual data is much larger. I am not using a specified coordinate system.
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
from osgeo import gdal
# X-coordinates (in meters)
dataX = np.random.uniform(-50,50,100)
# Y-coordinates (in meters)
dataY = np.random.uniform(-50,50,100)
# some other variable
dataZ = np.random.uniform(-50,50,100)
# cell dimensions
binWidth = 4
binLength = 20
# min/max coordinates
xMin = min(dataX)
xMax = max(dataX)
yMin = min(dataY)
yMax = max(dataY)
# create histogram to compare with raster
plt.hist2d(dataX, dataY, bins = [np.arange(xMin, xMax+binWidth, binWidth), np.arange(yMin, yMax+binLength, binLength)], weights = dataZ)
# binning
statistic,x_edge,y_edge,binnumber = stats.binned_statistic_2d(dataX, dataY, dataZ, 'mean',bins = [np.arange(xMin, xMax+binWidth, binWidth), np.arange(yMin, yMax+binLength, binLength)] )
# get number of rows and columns
nrows,ncols = np.shape(statistic)
# create the raster
geotransform=(xMin,20,0,yMax,0, -4)
output_raster = gdal.GetDriverByName('GTiff').Create('testRaster.tif',ncols, nrows, 1 ,gdal.GDT_Float32)
output_raster.SetGeoTransform(geotransform)
output_raster.GetRasterBand(1).WriteArray(statistic)
output_raster.FlushCache()
plt.show()
I have played around with the 'SetGeoTransform' parameters (I am assuming this is where the problem is) multiple times but none of my attempts have gotten the desired results. For example, switching the pixel sizes by setting geotransform=(xMin,4,0,yMax,0, -20)
, yields the following result:
What exactly am I doing wrong ?
geotransform=(xMin,4,0,yMax,0, -20)
, but then the raster appears stretched. I have added an additional image to show it looks if I swap the pixel sizes as you suggested.nrows,ncols = np.shape(statistic)
b) I want the pixels to have a length of 20 m (Y-distance)and a width of 4 m(X-distance), as shown in the histogram.