I have an xyz file in which I'd like to apply the interpolation method of kriging
to via python.
In order to do this, I have come across the module PyKrige
and duplicated this example perfectly on my own machine https://colab.research.google.com/drive/1mZDVJ-RioBc4QHHpe2HBNJrx8sJ4-1Lw?pli=1#scrollTo=z0XaQWJq_gJ7
This gave me hope that the module could be easily implemented on my own data - sadly not.
My xyz file looks like so, and there are 36,000 sets of xyz points (973KB)
626814.86 233370.97 14.56
626814.43 233371.28 14.55
626814.01 233371.56 14.54
626813.57 233371.85 14.53
626813.10 233372.15 14.51
626812.66 233372.44 14.49
626812.23 233372.72 14.46
626811.78 233373.01 14.43
I have then altered the script from the example above to suit my data, by reading them into their 3 separate arrays and what not
import numpy as np
import matplotlib.pyplot as plt
import pykrige.kriging_tools as kt
from pykrige.ok import OrdinaryKriging
from pylab import *
x1 = np.loadtxt("c:/test.xyz", usecols=0, dtype='float')
y1 = np.loadtxt("c:/test.xyz", usecols=1, dtype='float')
z1 = np.loadtxt("c:/test.xyz", usecols=2, dtype='float')
x = list(x1) #may not be a necessary step, but I added it so that after each x coordinate there was a comma like the examples array. Script wasn't working even before doing this
y = list(y1)
z = list(z1)
cax = plt.scatter(x, y, c=z)
cbar = plt.colorbar(cax, fraction=0.03)
plt.title('Test Points')
OK = OrdinaryKriging(
x,
y,
z,
variogram_model='gaussian',
verbose=True,
enable_plotting=False,
nlags=5,
)
OK.variogram_model_parameters
gridx = np.arange(626000, 629094, 5, dtype='float64')
gridy = np.arange(231500, 234500, 5, dtype='float64')
zstar, ss = OK.execute("grid", gridx, gridy)
print(zstar.shape)
print(ss.shape)
cax = plt.imshow(zstar, extent=(626000, 629094, 231500, 234500), origin='lower')
plt.scatter(x, y, c='z', marker='.')
cbar=plt.colorbar(cax)
plt.title('Test')
savefig( 'C:/users/public/test.png' )
print('image saved')
My issue is that it gets to the execution of the Ordinary Kriging and then complains
numpy.core._exceptions._ArrayMemoryError: Unable to allocate 10.2 GiB for an array with shape (36936, 36936) and data type float64
Understandably, I thought that the issue was the 36,000 points being too many. So, I decided to try a smaller snippet of the data. Now when I say smaller, I mean 20 points (1KB!!!). THE SAME ERROR!!!!
Can someone please tell me what I'm doing wrong? Realistically, I only want the interpolated PNG, not the points and scatterplot