1 文本格式
using System;
namespace Legalsoft.Truffer
{
/// <summary>
/// 径向基函数插值
/// Object for radial basis function interpolation using n points in dim
/// dimensions.Call constructor once, then interp as many times as desired.
/// </summary>
public class RBF_interp
{
//public delegateFunc func { get; set; } = null;
private int dim { get; set; }
private int n { get; set; }
private double[,] pts { get; set; }
private double[] vals { get; set; }
private double[] w;
public RBF_fn fn;
private bool norm { get; set; }
public RBF_interp(double[,] ptss, double[] valss, RBF_fn func, bool nrbf = false)
{
this.dim = ptss.GetLength(1);
this.n = ptss.GetLength(0);
this.pts = ptss;
this.vals = valss;
this.w = new double[n];
this.fn = func;
this.norm = nrbf;
double[,] rbf = new double[n, n];
double[] rhs = new double[n];
for (int i = 0; i < n; i++)
{
double sum = 0.0;
for (int j = 0; j < n; j++)
{
//sum += (rbf[i, j] = fn.rbf(rad(pts[i, 0], pts[j, 0])));
double r = Globals.dist(Globals.CopyFrom(i, pts), Globals.CopyFrom(j, pts));
rbf[i, j] = func.rbf(r);
sum += rbf[i, j];
}
if (norm)
{
rhs[i] = sum * vals[i];
}
else
{
rhs[i] = vals[i];
}
}
LUdcmp lu = new LUdcmp(rbf);
lu.solve( rhs, w);
}
public double interp(double[] pt)
{
double fval;
double sum = 0.0;
double sumw = 0.0;
if (pt.Length != dim)
{
throw new Exception("RBF_interp bad pt size");
}
for (int i = 0; i < n; i++)
{
//fval = fn.rbf(rad(pt[0], pts[i, 0]));
fval = fn.rbf(Globals.dist(pt, Globals.CopyFrom(i, pts)));
sumw += w[i] * fval;
sum += fval;
}
return norm ? sumw / sum : sumw;
}
/*
public double rad(double[] p1, double[] p2)
{
double sum = 0.0;
for (int i = 0; i < dim; i++)
{
sum += Globals.SQR(p1[i] - p2[i]);
}
if (sum <= float.Epsilon)
{
return 0.0;
}
return Math.Sqrt(sum);
}
*/
}
}
2 代码格式
using System;
namespace Legalsoft.Truffer
{
/// <summary>
/// 径向基函数插值
/// Object for radial basis function interpolation using n points in dim
/// dimensions.Call constructor once, then interp as many times as desired.
/// </summary>
public class RBF_interp
{
//public delegateFunc func { get; set; } = null;
private int dim { get; set; }
private int n { get; set; }
private double[,] pts { get; set; }
private double[] vals { get; set; }
private double[] w;
public RBF_fn fn;
private bool norm { get; set; }
public RBF_interp(double[,] ptss, double[] valss, RBF_fn func, bool nrbf = false)
{
this.dim = ptss.GetLength(1);
this.n = ptss.GetLength(0);
this.pts = ptss;
this.vals = valss;
this.w = new double[n];
this.fn = func;
this.norm = nrbf;
double[,] rbf = new double[n, n];
double[] rhs = new double[n];
for (int i = 0; i < n; i++)
{
double sum = 0.0;
for (int j = 0; j < n; j++)
{
//sum += (rbf[i, j] = fn.rbf(rad(pts[i, 0], pts[j, 0])));
double r = Globals.dist(Globals.CopyFrom(i, pts), Globals.CopyFrom(j, pts));
rbf[i, j] = func.rbf(r);
sum += rbf[i, j];
}
if (norm)
{
rhs[i] = sum * vals[i];
}
else
{
rhs[i] = vals[i];
}
}
LUdcmp lu = new LUdcmp(rbf);
lu.solve( rhs, w);
}
public double interp(double[] pt)
{
double fval;
double sum = 0.0;
double sumw = 0.0;
if (pt.Length != dim)
{
throw new Exception("RBF_interp bad pt size");
}
for (int i = 0; i < n; i++)
{
//fval = fn.rbf(rad(pt[0], pts[i, 0]));
fval = fn.rbf(Globals.dist(pt, Globals.CopyFrom(i, pts)));
sumw += w[i] * fval;
sum += fval;
}
return norm ? sumw / sum : sumw;
}
/*
public double rad(double[] p1, double[] p2)
{
double sum = 0.0;
for (int i = 0; i < dim; i++)
{
sum += Globals.SQR(p1[i] - p2[i]);
}
if (sum <= float.Epsilon)
{
return 0.0;
}
return Math.Sqrt(sum);
}
*/
}
}