Commit 3cadc40a authored by Nils Rummler's avatar Nils Rummler
Browse files

Final structur of directories.

parent 418b41be
......@@ -16,8 +16,8 @@
*
*/
#ifndef IW_ITA_BASE_SPLINE
#define IW_ITA_BASE_SPLINE
#ifndef IW_ITA_BASE_PIECEWISEPOLYNOMIAL
#define IW_ITA_BASE_PIECEWISEPOLYNOMIAL
#include <ITABaseDefinitions.h>
#include <vector>
......@@ -81,11 +81,6 @@ namespace ITABase
*/
std::vector<float> vdCoefficients;
};
//! Uses cubic splines to create a piecewise polynomial (order 3) for the given data pairs vdSupportingPoints and vdDataPoints. vdSupportingPoints need to be sorted in ascending order
ITA_BASE_API CPiecewisePolynomial CubicSpline(const std::vector<float>& vdSupportingPoints, const std::vector<float>& vdDataPoints);
//! Uses cubic splines to create a piecewise polynomial (order 3) for the given data pairs vdSupportingPoints and vdDataPoints. Evaluates that Polynom at vdCostumPoints. vdSupportingPoints need to be sorted in ascending order.
ITA_BASE_API std::vector<float> CubicSpline(const std::vector<float>& vdSupportingPoints, const std::vector<float>& vdDataPoints, const std::vector<float>& vdCostumPoints);
}
#endif // IW_ITA_BASE_SPLINE
#endif // IW_ITA_BASE_PIECEWISEPOLYNOMIAL
/*
* ----------------------------------------------------------------
*
* ITA core libs
* (c) Copyright Institute of Technical Acoustics (ITA)
* RWTH Aachen University, Germany, 2015-2020
*
* ----------------------------------------------------------------
* ____ __________ _______
* // / //__ ___/ // _ |
* // / // / // /_| |
* // / // / // ___ |
* //__/ //__/ //__/ |__|
*
* ----------------------------------------------------------------
*
*/
#ifndef IW_ITA_BASE_SPLINE
#define IW_ITA_BASE_SPLINE
#include <ITABaseDefinitions.h>
#include <ITABase/Math/PiecewisePolynomial.h>
#include <vector>
namespace ITABase
{
//! Uses cubic splines to create a piecewise polynomial (order 3) for the given data pairs vdSupportingPoints and vdDataPoints. vdSupportingPoints need to be sorted in ascending order
ITA_BASE_API CPiecewisePolynomial CubicSpline(const std::vector<float>& vdSupportingPoints, const std::vector<float>& vdDataPoints);
//! Uses cubic splines to create a piecewise polynomial (order 3) for the given data pairs vdSupportingPoints and vdDataPoints. Evaluates that Polynom at vdCostumPoints. vdSupportingPoints need to be sorted in ascending order.
ITA_BASE_API std::vector<float> CubicSpline(const std::vector<float>& vdSupportingPoints, const std::vector<float>& vdDataPoints, const std::vector<float>& vdCostumPoints);
}
#endif // IW_ITA_BASE_SPLINE
......@@ -6,6 +6,8 @@ set( SubDirs )
set( DirFiles
LinearAlgebra.h
PiecewisePolynomial.h
Spline.h
)
......
......@@ -6,7 +6,7 @@ set( SubDirs Math )
set( DirFiles
ITAProgress.h
Spline.h
)
if( ITA_BASE_WITH_JSON_SUPPORT )
......
#include <ITABase/Math/PiecewisePolynomial.h>
#include <ITAException.h>
#include <string>
#include <algorithm>
#include <cmath>
#include <ITABase/Math/LinearAlgebra.h>
using namespace ITABase;
CPiecewisePolynomial::CPiecewisePolynomial(const std::vector<float>& vdBreakPoints, const std::vector<float>& vdCoeffs, const int iOrder)
: iOrder(iOrder)
, vdBreakPoints(vdBreakPoints)
, vdCoefficients(vdCoeffs)
{
if (iOrder < 0)
ITA_EXCEPT_INVALID_PARAMETER("Polynom order must be greater or equal zero.");
const int iMinimumNBreakPoints = std::max(2, iOrder + 1); // At least one interval
if ( vdBreakPoints.size() < iMinimumNBreakPoints )
ITA_EXCEPT_INVALID_PARAMETER("Number of break points not sufficient for polynom order: Expected " + std::to_string(iMinimumNBreakPoints) + " but got " + std::to_string(vdBreakPoints.size()) + ".");
const int nCoeffs = NumIntervals() * (iOrder + 1);
if (vdCoeffs.size() != nCoeffs)
ITA_EXCEPT_INVALID_PARAMETER("Invalid number of coefficients: Expected " + std::to_string(nCoeffs) + " but got " + std::to_string(vdCoeffs.size()) + ".");
}
float CPiecewisePolynomial::Value(const float xValue) const
{
// Check if x-value is inside boundaries
// Find index of interval
int iInterval = IntervalIndex(xValue);
//TODO: get corresponding coefficients and calculate values
std::vector<float>::const_iterator itCoef = CoefficientsOfInterval(iInterval);
//f(x) = a1(x−x1)^3 + b1(x−x1)^2 + c1(x−x1) + d1
float x0 = vdBreakPoints[iInterval];
if (iOrder == 3)
return itCoef[0] * pow(xValue - x0, 3) + itCoef[1] * pow(xValue - x0, 2) + itCoef[2] * (xValue - x0) + itCoef[3];
else
ITA_EXCEPT_INVALID_PARAMETER("currently only splines of order 3");
}
std::vector<float> ITABase::CPiecewisePolynomial::Value(const std::vector<float>& vdValues) const
{
std::vector<float> result(vdValues.size());
for (int i = 0; i < vdValues.size(); i++)
result[i] = Value(vdValues[i]);
return result;
}
int ITABase::CPiecewisePolynomial::IntervalIndex(const float xValue) const
{
if (xValue < vdBreakPoints[0] || xValue > vdBreakPoints[NumIntervals()])
ITA_EXCEPT_INVALID_PARAMETER("xValue is out of bounds");
for (int i = 1; i < vdBreakPoints.size(); i++) {
if (vdBreakPoints[i] > xValue)
return i - 1;
}
}
CPiecewisePolynomial CPiecewisePolynomial::Derivation() const
{
//TODO: Check if this works correctly
const int iNewOrder = std::max(0, iOrder - 1); // Minimum order is zero
const int iNewNumCoefficients = std::max(1, NumCoefficients() - 1); // At least one coefficient
std::vector<float> vdNewCoeffs( NumIntervals() * iNewNumCoefficients, 0.0);
for (int idxPiece = 0; idxPiece < NumIntervals(); idxPiece++)
{
for (int idxCoeff = 0; idxCoeff < iOrder; idxCoeff++)
{
const int iCurrentOrder = iOrder - idxCoeff;
const int idxOld = idxPiece * NumCoefficients() + idxCoeff;
const int idxNew = idxPiece * iNewNumCoefficients + idxCoeff;
vdNewCoeffs[idxNew] = vdCoefficients[idxOld] * iCurrentOrder;
}
}
return CPiecewisePolynomial(vdBreakPoints, vdNewCoeffs, iNewOrder);
}
std::vector<float>::const_iterator ITABase::CPiecewisePolynomial::CoefficientsOfInterval(const int idxInterval) const
{
if (idxInterval < 0 || idxInterval >= NumIntervals())
ITA_EXCEPT_INVALID_PARAMETER("Index for interval is out of bounds.");
return vdCoefficients.begin() + idxInterval * NumCoefficients();
}
\ No newline at end of file
#include <ITABase/Spline.h>
#include <ITABase/Math/Spline.h>
#include <ITAException.h>
#include <string>
#include <algorithm>
#include <cmath>
#include <ITABase/Math/LinearAlgebra.h>
using namespace ITABase;
CPiecewisePolynomial::CPiecewisePolynomial(const std::vector<float>& vdBreakPoints, const std::vector<float>& vdCoeffs, const int iOrder)
: iOrder(iOrder)
, vdBreakPoints(vdBreakPoints)
, vdCoefficients(vdCoeffs)
{
if (iOrder < 0)
ITA_EXCEPT_INVALID_PARAMETER("Polynom order must be greater or equal zero.");
const int iMinimumNBreakPoints = std::max(2, iOrder + 1); // At least one interval
if ( vdBreakPoints.size() < iMinimumNBreakPoints )
ITA_EXCEPT_INVALID_PARAMETER("Number of break points not sufficient for polynom order: Expected " + std::to_string(iMinimumNBreakPoints) + " but got " + std::to_string(vdBreakPoints.size()) + ".");
const int nCoeffs = NumIntervals() * (iOrder + 1);
if (vdCoeffs.size() != nCoeffs)
ITA_EXCEPT_INVALID_PARAMETER("Invalid number of coefficients: Expected " + std::to_string(nCoeffs) + " but got " + std::to_string(vdCoeffs.size()) + ".");
}
float CPiecewisePolynomial::Value(const float xValue) const
{
// Check if x-value is inside boundaries
// Find index of interval
int iInterval = IntervalIndex(xValue);
//TODO: get corresponding coefficients and calculate values
std::vector<float>::const_iterator itCoef = CoefficientsOfInterval(iInterval);
//f(x) = a1(x−x1)^3 + b1(x−x1)^2 + c1(x−x1) + d1
float x0 = vdBreakPoints[iInterval];
if (iOrder == 3)
return itCoef[0] * pow(xValue - x0, 3) + itCoef[1] * pow(xValue - x0, 2) + itCoef[2] * (xValue - x0) + itCoef[3];
else
ITA_EXCEPT_INVALID_PARAMETER("currently only splines of order 3");
}
std::vector<float> ITABase::CPiecewisePolynomial::Value(const std::vector<float>& vdValues) const
{
std::vector<float> result(vdValues.size());
for (int i = 0; i < vdValues.size(); i++)
result[i] = Value(vdValues[i]);
return result;
}
int ITABase::CPiecewisePolynomial::IntervalIndex(const float xValue) const
{
if (xValue < vdBreakPoints[0] || xValue > vdBreakPoints[NumIntervals()])
ITA_EXCEPT_INVALID_PARAMETER("xValue is out of bounds");
for (int i = 1; i < vdBreakPoints.size(); i++) {
if (vdBreakPoints[i] > xValue)
return i - 1;
}
}
CPiecewisePolynomial CPiecewisePolynomial::Derivation() const
{
//TODO: Check if this works correctly
const int iNewOrder = std::max(0, iOrder - 1); // Minimum order is zero
const int iNewNumCoefficients = std::max(1, NumCoefficients() - 1); // At least one coefficient
std::vector<float> vdNewCoeffs( NumIntervals() * iNewNumCoefficients, 0.0);
for (int idxPiece = 0; idxPiece < NumIntervals(); idxPiece++)
{
for (int idxCoeff = 0; idxCoeff < iOrder; idxCoeff++)
{
const int iCurrentOrder = iOrder - idxCoeff;
const int idxOld = idxPiece * NumCoefficients() + idxCoeff;
const int idxNew = idxPiece * iNewNumCoefficients + idxCoeff;
vdNewCoeffs[idxNew] = vdCoefficients[idxOld] * iCurrentOrder;
}
}
return CPiecewisePolynomial(vdBreakPoints, vdNewCoeffs, iNewOrder);
}
std::vector<float>::const_iterator ITABase::CPiecewisePolynomial::CoefficientsOfInterval(const int idxInterval) const
{
if (idxInterval < 0 || idxInterval >= NumIntervals())
ITA_EXCEPT_INVALID_PARAMETER("Index for interval is out of bounds.");
return vdCoefficients.begin() + idxInterval * NumCoefficients();
}
ITA_BASE_API CPiecewisePolynomial ITABase::CubicSpline(const std::vector<float>& vdSupportingPoints, const std::vector<float>& vdDataPoints)
{
/*setUp the bandmatrix M with the the h values. The vector with the coefficent b is called b.
......
......@@ -6,6 +6,8 @@ set( SubDirs )
set( DirFiles
LinearAlgebra.cpp
PiecewisePolynomial.cpp
Spline.cpp
)
......
......@@ -6,7 +6,6 @@ set( SubDirs Math )
set( DirFiles
ITAProgress.cpp
Spline.cpp
)
if( ITA_BASE_WITH_JSON_SUPPORT )
......
......@@ -19,8 +19,8 @@
*/
#include <ITABase/Spline.h>
#include <ITABase/Math/Spline.h>
#include <ITAStopWatch.h>
#include <iostream>
//using namespace std;
......@@ -44,10 +44,10 @@ int main(int iNumInArgs, char* pcInArgs[])
int iPolynomialOrder = 3;
IW_ITA_BASE_SPLINE::ITABase::CPiecewisePolynomial myPoly1 = IW_ITA_BASE_SPLINE::ITABase::Spline(vdSupportingPoints,vdDataPoints,iPolynomialOrder);
if MyCompare(std::cout << "Quick and dirty test: "<< myPoly1.Value(5)<<" -hier soll grob 10 rauskommen;)";
*/
float epsilon = 0.1;
//Matlab spline. DataPoints sind aus randi
std::vector<float> vdSupportingPoints{ -10, - 9, - 8, - 7, - 6, - 5, - 4, - 3, - 2, - 1, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 };
std::vector<float> vdDataPoints = { 13,14,2,14,10,2,5,9,15,15,3,15,11,4,0,-1,-5,-15,-11,-3,2 };
std::vector<float> vdX = { -10, float(1.5), float(-7.898),3, float(3.3), float(9.9),-2 };
std::vector<float> vdResult = { float(13.0000), float(17.4325), float(2.4282), float(8.0000), float(9.4143), float(11.4819), float(15.0000) };
......@@ -56,17 +56,31 @@ int main(int iNumInArgs, char* pcInArgs[])
std::cout << "MATLAB SPLINE TEST 1" << std::endl;
if (MyCompare(vdResult, myPoly1.Value(vdX),epsilon))
std::cout << "Test 1 bestanden" << std::endl;
*/
//Test 2
std::vector<float> vdResult2 = { float(77.0000), float(7.9485), float(18.7512), float(50.0000), float(70.7647), float(36.3051), float(28.0000)
*///Test 2
std::vector<float> vdSupportingPoints{ -10, -9, -8, -7, -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 };
std::vector<float> vdDataPoints2{ 77, 80, 19, 49, 45, 65, 71, 76, 28, 68, 66, 17,12, 50, 96, 35, 59, 23, 76, 26, 51 };
std::vector<float> vdResult2 = { float(77.0000), float(7.9485), float(18.7512), float(50.0000), float(70.7647), float(36.3051), float(28.0000) };
iPolynomialOrder = 3;
IW_ITA_BASE_SPLINE::ITABase::CPiecewisePolynomial myPoly2 = IW_ITA_BASE_SPLINE::ITABase::Spline(vdSupportingPoints, vdDataPoints2, iPolynomialOrder);
std::cout << "MATLAB SPLINE TEST 2" << std::endl;
if (MyCompare(vdResult2, myPoly2.Value(vdX),epsilon))
std::cout << "Test 2 bestanden" << std::endl;
std::vector<float> vdX = { -10, float(1.5), float(-7.898),float(3.3), float(3.3), float(9.9),float(-2.5),float(8.9),float(0.01),float(-4.3) };
IW_ITA_BASE_SPLINE::ITABase::CPiecewisePolynomial myPoly= ITABase::CubicSpline(vdSupportingPoints, vdDataPoints2);
std::vector<float> result(10);
ITAStopWatch myWatch = ITAStopWatch();
std::cout << "Rueckgabe von piecwisePolynomial Objekten" << std::endl;
for (int k = 0; k < 100; k++) {
myWatch.start();
myPoly = ITABase::CubicSpline(vdSupportingPoints, vdDataPoints2);
myWatch.stop();
}
std::cout << myWatch << std::endl;
myWatch.reset();
std::cout << "Rueckgabe von vector Objekten" << std::endl;
for (int k = 0; k < 100; k++) {
myWatch.start();
result = ITABase::CubicSpline(vdSupportingPoints, vdDataPoints2, vdX);
myWatch.stop();
}
std::cout << myWatch << std::endl;
return 0;
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment