kriging克⾥⾦插值以及前端渲染jS代码部分解释
.html中的⽅法调⽤
//训练使⽤⾼斯过程与贝叶斯先验
let ain(positions.map(pos=>pos[2]),positions.map(pos=>pos[0]),positions.map(pos=>pos[1]),params.krigingModel,params.krigingSig ma2,params.krigingAlpha);
//⽹格矩阵或轮廓路径
let id(polygons,variogram,(extent[2]-extent[0])/1000);
//在DOM上绘图.
//Canvas是HTML5提供的⼀个标签,我们可以在这个盒⼦区域绘画
kriging.plot(canvas,grid,[extent[0],extent[2]],[extent[1],extent[3]],colors);
kriging-original.js⽂件中的部分代码解释
// Extend the Array class
// Array.prototype.max重写数组原型链
//表⽰取得最⼤值
Array.prototype.max =function(){
//apply()⽅法接受的是⼀个参数数组
//返回⼀个最⼤值的数组
return Math.max.apply(null,this);
};
//这⾥表⽰取得最⼩值
Array.prototype.min =function(){
//返回⼀个最⼩值
return Math.min.apply(null,this);
};
//这⾥表⽰算平均数
an =function(){
var i, sum;
for(i =0, sum =0; i <this.length; i++)
sum +=this[i];
return sum /this.length;
};
p =function(n){
//返回⼀个长度为n的数组,且每⼀个元素都被赋值成undefined
return Array.apply(null,new Array(n))
.map(Number.prototype.valueOf,this[0]);
//Number.prototype.valueOf()⽅法返回数值对象的原始值
};
Array.prototype.pip =function(x, y){
js购物车结算代码var i, j, c =false;
for(i =0, j =this.length -1; i <this.length; j = i++){
if(((this[i][1]> y)!=(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;
}
var kriging =function(){
var kriging ={};
// Matrix algebra矩阵代数
kriging_matrix_diag =function(c, n){
//新建⼀个n*n的矩阵
var i, Z =[0].rep(n * n);
//循环赋值c给Z矩阵的每⼀元素
for(i =0; i < n; i++) Z[i * n + i]= c;
return Z;
};
};
//将这个矩阵变为转置阵,也就是将元素颠倒顺序
kriging_matrix_transpose =function(X, n, m){
var i, j, Z =Array(m * n);
for(i =0; i < n; i++)
for(j =0; j < m; j++)
Z[j * n + i]= X[i * m + j];
return Z;
};
/
/再次改变数值,把c给每⼀个⼆维元素赋值
kriging_matrix_scale =function(X, c, n, m){
var i, j;
for(i =0; i < n; i++)
for(j =0; j < m; j++)
X[i * m + j]*= c;
};
//添加的⽅法
kriging_matrix_add =function(X, Y, n, m){
//新建⼀个m*n的矩阵Z
var i, j, Z =Array(n * m);
for(i =0; i < n; i++)
for(j =0; j < m; j++)
//将X和Y矩阵相加合并成⼀个矩阵
Z[i * m + j]= X[i * m + j]+ Y[i * m + j];
//返回⼀个Z矩阵
return Z;
};
// Naive matrix multiplication
//简单的矩阵乘法,矩阵和矩阵的乘法
//也就是前⼀个矩阵中的⾏乘以后⼀个矩阵中的列
kriging_matrix_multiply =function(X, Y, n, m, p){
var i, j, k, Z =Array(n * p);
for(i =0; i < n; i++){
for(j =0; j < p; j++){
Z[i * p + j]=0;
for(k =0; k < m; k++)
Z[i * p + j]+= X[i * m + k]* Y[k * p + j];
}
}
return Z;
};
// Cholesky decomposition
/
/柯列斯基分解,这是⼀种将正定矩阵分解为上三⾓矩阵和下三⾓矩阵的⽅法,//在优化矩阵计算的时候会⽤到的⼀种技巧
//也就是,下⾯左边为下三⾓,右边为上三⾓
//100000 123456
//120000 023456
//123000 003456
//123400 000456
//123450 000056
//123456 000006
kriging_matrix_chol =function(X, n){
var i, j, k, sum, p =Array(n);
for(i =0; i < n; i++) p[i]= X[i * n + i];
for(i =0; i < n; i++){
for(j =0; j < i; j++)
p[i]-= X[i * n + j]* X[i * n + j];
if(p[i]<=0)return false;
p[i]= Math.sqrt(p[i]);
for(j = i +1; j < n; j++){
for(k =0; k < i; k++)
X[j * n + i]-= X[j * n + k]* X[i * n + k];
X[j * n + i]/= p[i];
}
}
for(i =0; i < n; i++) X[i * n + i]= p[i];
for(i =0; i < n; i++) X[i * n + i]= p[i];
return true;
};
// Inversion of cholesky decomposition
//⽤斯基分解求矩阵的逆
kriging_matrix_chol2inv =function(X, n){
var i, j, k, sum;
for(i =0; i < n; i++){
X[i * n + i]=1/ X[i * n + i];
for(j = i +1; j < n; j++){
sum =0;
for(k = i; k < j; k++)
sum -= X[j * n + k]* X[k * n + i];
X[j * n + i]= sum / X[j * n + j];
}
}
for(i =0; i < n; i++)
for(j = i +1; j < n; j++)
X[i * n + j]=0;
for(i =0; i < n; i++){
X[i * n + i]*= X[i * n + i];
for(k = i +1; k < n; k++)
X[i * n + i]+= X[k * n + i]* X[k * n + i];
for(j = i +1; j < n; j++)
for(k = j; k < n; k++)
X[i * n + j]+= X[k * n + i]* X[k * n + j];
}
for(i =0; i < n; i++)
for(j =0; j < i; j++)
X[i * n + j]= X[j * n + i];
};
// Inversion via gauss-jordan elimination
/
/⽤⾼斯-约当消去法求逆,它的速度不是最快的,但是它⾮常稳定
//如果A是求解矩阵,那么求A的逆矩阵则为
//⽤A矩阵右边乘以单位矩阵I(与A同⾏同列值为1的单位矩阵)
//公式为A*I=I*B,(等号右边要同时变化),也就是⼀个矩阵右边乘以单位矩阵化为,//左边单位矩阵乘以B,则B就是A矩阵的逆
kriging_matrix_solve =function(X, n){
var m = n;
var b =Array(n * n);
var indxc =Array(n);
var indxr =Array(n);
var ipiv =Array(n);
var i, icol, irow, j, k, l, ll;
var big, dum, pivinv, temp;
for(i =0; i < n; i++)
for(j =0; j < n; j++){
if(i == j) b[i * n + j]=1;
else b[i * n + j]=0;
}
for(j =0; j < n; j++) ipiv[j]=0;
for(i =0; i < n; i++){
big =0;
for(j =0; j < n; j++){
if(ipiv[j]!=1){
for(k =0; k < n; k++){
if(ipiv[k]==0){
if(Math.abs(X[j * n + k])>= big){
big = Math.abs(X[j * n + k]);
irow = j;
icol = k;
}
}
}
}
}
}
++(ipiv[icol]);
if(irow != icol){
for(l =0; l < n; l++){
temp = X[irow * n + l];
X[irow * n + l]= X[icol * n + l];
X[icol * n + l]= temp;
}
for(l =0; l < m; l++){
temp = b[irow * n + l];
b[irow * n + l]= b[icol * n + l];
b[icol * n + l]= temp;
}
}
indxr[i]= irow;
indxc[i]= icol;
if(X[icol * n + icol]==0)return false;// Singular
pivinv =1/ X[icol * n + icol];
X[icol * n + icol]=1;
for(l =0; l < n; l++) X[icol * n + l]*= pivinv;
for(l =0; l < m; l++) b[icol * n + l]*= pivinv;
for(ll =0; ll < n; ll++){
if(ll != icol){
dum = X[ll * n + icol];
X[ll * n + icol]=0;
for(l =0; l < n; l++) X[ll * n + l]-= X[icol * n + l]* dum;
for(l =0; l < m; l++) b[ll * n + l]-= b[icol * n + l]* dum;
}
}
}
for(l =(n -1); l >=0; l--)
if(indxr[l]!= indxc[l]){
for(k =0; k < n; k++){
temp = X[k * n + indxr[l]];
X[k * n + indxr[l]]= X[k * n + indxc[l]];
X[k * n + indxc[l]]= temp;
}
}
return true;
}
// Variogram models
//变差函数模型
//变差函数⾼斯
kriging_variogram_gaussian =function(h, nugget, range, sill, A){ return nugget +((sill - nugget)/ range)*
(1.p(-(1.0/ A)* Math.pow(h / range,2)));
};
//变差函数指数
kriging_variogram_exponential =function(h, nugget, range, sill, A){ return nugget +((sill - nugget)/ range)*
(1.p(-(1.0/ A)*(h / range)));
};
//变差函数的球形
kriging_variogram_spherical =function(h, nugget, range, sill, A){ if(h > range)return nugget +(sill - nugget)/ range;
return nugget +((sill - nugget)/ range)*
(1.5*(h / range)-0.5* Math.pow(h / range,3));
(1.5*(h / range)-0.5* Math.pow(h / range,3));
};
// Train using gaussian processes with bayesian priors
//训练使⽤⾼斯过程与贝叶斯先验
//ain(t, x, y, model, sigma2, alpha):
//使⽤gaussian、exponential或spherical模型对数据集进⾏训练,返回的是⼀个variogram对象; ain =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":
break;
case"exponential":
break;
case"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; i < n; i++)
for(j =0; j < i; j++, k++){
distance[k]=Array(2);
distance[k][0]= Math.pow(
Math.pow(x[i]- x[j],2)+
Math.pow(y[i]- y[j],2),0.5);
distance[k][1]= Math.abs(t[i]- t[j]);
}
distance.sort(function(a, b){return a[0]- b[0];});
variogram.range = distance[(n * n - n)/2-1][0];
// Bin lag distance
//本滞后距离
var lags =((n * n - n)/2)>30?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 < lags; l++){
lag[l]= distance[l][0];
semi[l]= distance[l][1];
}
}else{
for(i =0, j =0, k =0, l =0; i < lags && j <((n * n - n)/2); i++, k =0){
while(distance[j][0]<=((i +1)* tolerance)){
lag[l]+= distance[j][0];
semi[l]+= distance[j][1];
j++;
k++;
if(j >=((n * n - n)/2))break;
}
if(k >0){
版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。
发表评论