kriging = { // Extend the Array class Array.prototype.max = function() { return Math.max.apply(null, this); }; Array.prototype.min = function() { return Math.min.apply(null, this); }; Array.prototype.mean = function() { var i, sum; for(i=0,sum=0;iy) != (this[j][1]>y)) && (x<(this[j][0]-this[i][0]) * (y-this[i][1]) / (this[j][1]-this[i][1]) + this[i][0]) ) { c = !c; } } return c; } const kriging = {}; // Matrix algebra const kriging_matrix_diag = function(c, n) { var i, Z = [0].rep(n*n); for(i=0;i=big) { big = Math.abs(X[j*n+k]); irow = j; icol = k; } } } } } ++(ipiv[icol]); if(irow!=icol) { for(l=0;l=0;l--) if(indxr[l]!=indxc[l]) { for(k=0;krange) return nugget + (sill-nugget)/range; return nugget + ((sill-nugget)/range)* ( 1.5*(h/range) - 0.5*Math.pow(h/range, 3) ); }; // Train using gaussian processes with bayesian priors kriging.train = function(t, x, y, model, sigma2, alpha) { var variogram = { t : t, x : x, y : y, nugget : 0.0, range : 0.0, sill : 0.0, A : 1/3, n : 0 }; switch(model) { case "gaussian": variogram.model = kriging_variogram_gaussian; break; case "exponential": variogram.model = kriging_variogram_exponential; break; case "spherical": variogram.model = kriging_variogram_spherical; break; }; // Lag distance/semivariance var i, j, k, l, n = t.length; var distance = Array((n*n-n)/2); for(i=0,k=0;i30?30:(n*n-n)/2; var tolerance = variogram.range/lags; var lag = [0].rep(lags); var semi = [0].rep(lags); if(lags<30) { for(l=0;l=((n*n-n)/2)) break; } if(k>0) { lag[l] /= k; semi[l] /= k; l++; } } if(l<2) return variogram; // Error: Not enough points } // Feature transformation n = l; variogram.range = lag[n-1]-lag[0]; var X = [1].rep(2*n); var Y = Array(n); var A = variogram.A; for(i=0;ixlim[1]) xlim[1] = polygons[i][j][0]; if(polygons[i][j][1]ylim[1]) ylim[1] = polygons[i][j][1]; } // Alloc for O(n^2) space var xtarget, ytarget; var a = Array(2), b = Array(2); var lxlim = Array(2); // Local dimensions var lylim = Array(2); // Local dimensions var x = Math.ceil((xlim[1]-xlim[0])/width); var y = Math.ceil((ylim[1]-ylim[0])/width); var A = Array(x+1); for(i=0;i<=x;i++) A[i] = Array(y+1); for(i=0;ilxlim[1]) lxlim[1] = polygons[i][j][0]; if(polygons[i][j][1]lylim[1]) lylim[1] = polygons[i][j][1]; } // Loop through polygon subspace a[0] = Math.floor(((lxlim[0]-((lxlim[0]-xlim[0])%width)) - xlim[0])/width); a[1] = Math.ceil(((lxlim[1]-((lxlim[1]-xlim[1])%width)) - xlim[0])/width); b[0] = Math.floor(((lylim[0]-((lylim[0]-ylim[0])%width)) - ylim[0])/width); b[1] = Math.ceil(((lylim[1]-((lylim[1]-ylim[1])%width)) - ylim[0])/width); for(j=a[0];j<=a[1];j++) for(k=b[0];k<=b[1];k++) { xtarget = xlim[0] + j*width; ytarget = ylim[0] + k*width; if(polygons[i].pip(xtarget, ytarget)) A[j][k] = kriging.predict(xtarget, ytarget, variogram); } } A.xlim = xlim; A.ylim = ylim; A.zlim = [variogram.t.min(), variogram.t.max()]; A.width = width; return A; }; kriging.contour = function(value, polygons, variogram) { }; // Plotting on the DOM kriging.plot = function(canvas, grid, xlim, ylim, colors) { // Clear screen var ctx = canvas.getContext("2d"); ctx.clearRect(0, 0, canvas.width, canvas.height); // Starting boundaries var range = [xlim[1]-xlim[0], ylim[1]-ylim[0], grid.zlim[1]-grid.zlim[0]]; var i, j, x, y, z; var n = grid.length; var m = grid[0].length; var wx = Math.ceil(grid.width*canvas.width/(xlim[1]-xlim[0])); var wy = Math.ceil(grid.width*canvas.height/(ylim[1]-ylim[0])); for(i=0;i1.0) z = 1.0; ctx.fillStyle = colors[Math.floor((colors.length-1)*z)]; ctx.fillRect(Math.round(x-wx/2), Math.round(y-wy/2), wx, wy); } }; return kriging; }