App
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
app::CubicSpline Class Reference

#include <CubicSpline.hpp>

Public Member Functions

 CubicSpline (double const *x, double const *y, unsigned count, double y1_0=1.e99, double y1_n=1.e99)
 
int isOK () const
 
double getValueAt (double x)
 
virtual ~CubicSpline ()
 

Protected Member Functions

int locate (const double x)
 
int hunt (const double x)
 
double interpolate (int jlo, double x) const
 
int y2ValsInit (const double *xv, const double *yv, double y1_0, double y1_n)
 

Protected Attributes

unsigned valsCount
 velikost dat x
 
unsigned lastVal
 
int jsav
 ??? Poslední nalezená souřadnice
 
int cor
 ??? bool - jestli je |jl - jsav| > dj;, udává, jestli se aktivuje hunt() nebo locate()
 
int dj
 ??? Celkov2 heuristika pro vyhledávání koeficientů
 
int ascnd
 ??? Zbytečně se počítalo v locate...
 
double * xVals
 
double * yVals
 xx, yy
 
double * y2Vals
 
int initOK
 

Detailed Description

Definition at line 14 of file CubicSpline.hpp.

Constructor & Destructor Documentation

app::CubicSpline::CubicSpline ( double const *  x,
double const *  y,
unsigned  count,
double  y1_0 = 1.e99,
double  y1_n = 1.e99 
)

Definition at line 6 of file CubicSpline.cpp.

{
if (count < 2) { initOK = 0; return; }
valsCount = count;
lastVal = count - 1;
jsav = 0; cor = 0;
if (x == NULL) { initOK = 0; return; }
if (y == NULL) { initOK = 0; return; }
xVals = (double *) malloc((unsigned long int)(valsCount * sizeof(double))); // ANDROID_SPECIFIC přetypování v malloc
if (xVals == NULL) { initOK = 0; return; }
memcpy(xVals, x, (unsigned long int)(count * sizeof(double))); // xVals = x;
yVals = (double *) malloc((unsigned long int)(valsCount * sizeof(double))); // ANDROID_SPECIFIC přetypování v malloc
if (yVals == NULL) { initOK = 0; free(xVals); return; }
memcpy(yVals, y, (unsigned long int)(count * sizeof(double))); // yVals = y; // ANDROID_SPECIFIC přetypování v malloc
y2Vals = (double *) malloc((unsigned long int)(valsCount * sizeof(double))); // ANDROID_SPECIFIC přetypování v malloc
if (y2Vals == NULL) { initOK = 0; free(xVals); free(yVals); return; }
ascnd = (xVals[lastVal] >= xVals[0]);
dj = pow(count, 0.25);
dj = (dj > 1) ? 1 : dj; //MIN(1,(int)pow((Doub)n,0.25)); // WTF???
if ( y2ValsInit(xVals, yVals, y1_0, y1_n) ) { initOK = 0; free(xVals); free(yVals); free(y2Vals); return; }
initOK = 1;
}
app::CubicSpline::~CubicSpline ( )
virtual

Definition at line 33 of file CubicSpline.cpp.

{
if (xVals != NULL) { free(xVals); }
if (yVals != NULL) { free(yVals); }
if (y2Vals != NULL) { free(y2Vals); }
}

Member Function Documentation

double app::CubicSpline::getValueAt ( double  x)

Definition at line 39 of file CubicSpline.cpp.

{
if (initOK) { // Při inicializaci nenastal žádný problém. Díky tomuto zmizí testování v ostatních metodách
if (x <= xVals[0]) {
return yVals[0];
} else if (x >= xVals[lastVal]) {
return yVals[lastVal];
}
int jlo = cor ? hunt(x) : locate(x); // Prom je zbytečná, lze rovnou volat s tímto funkci
/*
if (jlo < 0) {
return 0; // Proč? Co to znamená? NESMYSL - nemůže vrátit menší, zajistěno při návratu
}
*/
return interpolate(jlo, x);
} else {
return 0;
}
}
int app::CubicSpline::hunt ( const double  x)
protected

Definition at line 82 of file CubicSpline.cpp.

{
/* NELZE
if (valsCount < 2) {
return -1; // throw("hunt size error");
}
*/
int jl = jsav;
int jm, ju;
int inc = 1;
if (jl < 0 || jl > lastVal) {
jl = 0;
ju = lastVal;
} else {
if (x >= xVals[jl] == ascnd) {
while (1) {
ju = jl + inc;
if (ju >= lastVal) {
ju = lastVal;
break;
} else if (x < xVals[ju] == ascnd)
break;
else {
jl = ju;
inc += inc;
}
}
} else {
ju = jl;
while (1) {
jl = jl - inc;
if (jl <= 0) {
jl = 0;
break;
} else if (x >= xVals[jl] == ascnd)
break;
else {
ju = jl;
inc += inc;
}
}
}
}
while (ju - jl > 1) {
jm = (ju + jl) >> 1;
if (x >= xVals[jm] == ascnd)
jl = jm;
else
ju = jm;
}
cor = abs(jl - jsav) > dj ? 0 : 1;
jsav = jl;
int retVal = valsCount - 2;
retVal = (retVal < jl) ? retVal : jl;
return (retVal > 0) ? retVal : 0; // MAX(0,MIN(valsCount-2,jl)); // MAX(0,MIN(n-mm,jl-((mm-2)>>1))); // n = count, mm = kernelWidth
}
double app::CubicSpline::interpolate ( int  jlo,
double  x 
) const
protected

Definition at line 174 of file CubicSpline.cpp.

{
int klo = jl, khi = jl + 1;
double y, h, b, a;
h = xVals[khi] - xVals[klo];
// if (h == 0.0) { return 0; } //throw("Bad input to routine splint"); // Toto by nikdy nemelo nastat, protoze mam vzdy min. 2 body a vzdy vracim jl o 2 mensi, nez je max...
a = (xVals[khi] - x) / h;
b = (x - xVals[klo]) / h;
y = a * yVals[klo] + b * yVals[khi]
+ ((a * a * a - a) * y2Vals[klo] + (b * b * b - b) * y2Vals[khi]) * (h * h)
/ 6.0;
return y;
}
int app::CubicSpline::isOK ( ) const
inline

Definition at line 35 of file CubicSpline.hpp.

{ return initOK; };
int app::CubicSpline::locate ( const double  x)
protected

Definition at line 58 of file CubicSpline.cpp.

{
/* NELZE - pojistka v initOK
if (valsCount < 2) {
return -1; // throw("locate size error");
}
*/
int jl = 0;
int ju = lastVal; //
int jm = 0;
while (ju - jl > 1) {
jm = (ju + jl) >> 1;
if (x >= xVals[jm] == ascnd) {
jl = jm;
} else {
ju = jm;
}
}
cor = abs(jl - jsav) > dj ? 0 : 1;
jsav = jl;
int retVal = valsCount - 2;
retVal = (retVal < jl) ? retVal : jl;
return (retVal > 0) ? retVal : 0; // je kontrola na zápornost nutná? Asi ne.. Ale okrajová optimalizace
}
int app::CubicSpline::y2ValsInit ( const double *  xv,
const double *  yv,
double  y1_0,
double  y1_n 
)
protected

Definition at line 138 of file CubicSpline.cpp.

{
// Jako optimalizace doplnena lastVal, pro pochopeni viz kniha.
int i, k;
double p, qn, sig, un;
double * u = (double *) malloc ((unsigned long int)((lastVal) * sizeof(double))); // ANDROID_SPECIFIC přetypování v malloc
if (u == NULL) { return 1; }
if (y1_0 > 0.99e99) { // Natural spline
y2Vals[0] = u[0] = 0.0;
} else {
y2Vals[0] = -0.5;
u[0] = (3.0 / (xv[1] - xv[0]))
* ((yv[1] - yv[0]) / (xv[1] - xv[0]) - y1_0);
}
for (i = 1; i < lastVal; i++) {
sig = (xv[i] - xv[i - 1]) / (xv[i + 1] - xv[i - 1]);
p = sig * y2Vals[i - 1] + 2.0;
y2Vals[i] = (sig - 1.0) / p;
u[i] = (yv[i + 1] - yv[i]) / (xv[i + 1] - xv[i])
- (yv[i] - yv[i - 1]) / (xv[i] - xv[i - 1]);
u[i] = (6.0 * u[i] / (xv[i + 1] - xv[i - 1]) - sig * u[i - 1]) / p;
}
if (y1_n > 0.99e99) { // Natural spline
qn = un = 0.0;
} else {
qn = 0.5;
un = (3.0 / (xv[lastVal] - xv[valsCount - 2]))
* (y1_n - (yv[lastVal] - yv[valsCount - 2]) / (xv[lastVal] - xv[valsCount - 2]));
}
y2Vals[lastVal] = (un - qn * u[valsCount - 2]) / (qn * y2Vals[valsCount - 2] + 1.0);
for (k = valsCount - 2; k >= 0; k--) {
y2Vals[k] = y2Vals[k] * y2Vals[k + 1] + u[k];
}
free(u);
return 0;
}

Member Data Documentation

int app::CubicSpline::ascnd
protected

??? Zbytečně se počítalo v locate...

Definition at line 21 of file CubicSpline.hpp.

int app::CubicSpline::cor
protected

??? bool - jestli je |jl - jsav| > dj;, udává, jestli se aktivuje hunt() nebo locate()

Definition at line 19 of file CubicSpline.hpp.

int app::CubicSpline::dj
protected

??? Celkov2 heuristika pro vyhledávání koeficientů

Definition at line 20 of file CubicSpline.hpp.

int app::CubicSpline::initOK
protected

Definition at line 25 of file CubicSpline.hpp.

int app::CubicSpline::jsav
protected

??? Poslední nalezená souřadnice

Definition at line 18 of file CubicSpline.hpp.

unsigned app::CubicSpline::lastVal
protected

Definition at line 17 of file CubicSpline.hpp.

unsigned app::CubicSpline::valsCount
protected

velikost dat x

Definition at line 16 of file CubicSpline.hpp.

double* app::CubicSpline::xVals
protected

Definition at line 22 of file CubicSpline.hpp.

double* app::CubicSpline::y2Vals
protected

Definition at line 23 of file CubicSpline.hpp.

double * app::CubicSpline::yVals
protected

xx, yy

Definition at line 22 of file CubicSpline.hpp.


The documentation for this class was generated from the following files: