Skip to content
Snippets Groups Projects
Commit 45771093 authored by Carl Philipp Klemm's avatar Carl Philipp Klemm
Browse files

seperate r-series out in libdrt to make the api more self explanitory

parent 6cf3485c
No related branches found
No related tags found
No related merge requests found
...@@ -168,7 +168,8 @@ static Eigen::Matrix<fv, Eigen::Dynamic, 2> calcBounds(Eigen::VectorX<std::compl ...@@ -168,7 +168,8 @@ static Eigen::Matrix<fv, Eigen::Dynamic, 2> calcBounds(Eigen::VectorX<std::compl
} }
template<typename fv> template<typename fv>
Eigen::VectorX<fv> calcDrt(Eigen::VectorX<std::complex<fv>>& impedanceSpectra, Eigen::VectorX<fv>& omegaTensor, FitMetics& fm, const FitParameters& fp) Eigen::VectorX<fv> calcDrt(Eigen::VectorX<std::complex<fv>>& impedanceSpectra, Eigen::VectorX<fv>& omegaTensor,
FitMetics& fm, const FitParameters& fp, fv* rSeries)
{ {
Eigen::Matrix<fv, Eigen::Dynamic, Eigen::Dynamic> aMatrixImag = aImag<fv>(omegaTensor); Eigen::Matrix<fv, Eigen::Dynamic, Eigen::Dynamic> aMatrixImag = aImag<fv>(omegaTensor);
Eigen::Matrix<fv, Eigen::Dynamic, Eigen::Dynamic> aMatrixReal = aReal<fv>(omegaTensor); Eigen::Matrix<fv, Eigen::Dynamic, Eigen::Dynamic> aMatrixReal = aReal<fv>(omegaTensor);
...@@ -201,14 +202,17 @@ Eigen::VectorX<fv> calcDrt(Eigen::VectorX<std::complex<fv>>& impedanceSpectra, E ...@@ -201,14 +202,17 @@ Eigen::VectorX<fv> calcDrt(Eigen::VectorX<std::complex<fv>>& impedanceSpectra, E
throw drt_errror(std::string(ex.what())); throw drt_errror(std::string(ex.what()));
} }
return x; if(rSeries)
*rSeries = x[x.size()-1];
return x(Eigen::seq(0, x.size()-2));
} }
template Eigen::VectorX<double> calcDrt<double>(Eigen::VectorX<std::complex<double>>&, template Eigen::VectorX<double> calcDrt<double>(Eigen::VectorX<std::complex<double>>&,
Eigen::VectorX<double>&, FitMetics& fm, const FitParameters& fp); Eigen::VectorX<double>&, FitMetics& fm, const FitParameters& fp, double* rSeries);
template Eigen::VectorX<float> calcDrt<float>(Eigen::VectorX<std::complex<float>>&, template Eigen::VectorX<float> calcDrt<float>(Eigen::VectorX<std::complex<float>>&,
Eigen::VectorX<float>&, FitMetics& fm, const FitParameters& fp); Eigen::VectorX<float>&, FitMetics& fm, const FitParameters& fp, float* rSeries);
template<typename fv> template<typename fv>
fv testFn() fv testFn()
...@@ -220,21 +224,23 @@ fv testFn() ...@@ -220,21 +224,23 @@ fv testFn()
template double testFn<double>(); template double testFn<double>();
#ifdef USE_EISGEN #ifdef USE_EISGEN
std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, const std::vector<fvalue>& omegaVector, FitMetics& fm, const FitParameters& fp) std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, const std::vector<fvalue>& omegaVector,
FitMetics& fm, const FitParameters& fp, fvalue* rSeries)
{ {
Eigen::VectorX<std::complex<fvalue>> impedanceSpectra = eistoeigen(data); Eigen::VectorX<std::complex<fvalue>> impedanceSpectra = eistoeigen(data);
Eigen::VectorX<fvalue> omega = Eigen::VectorX<fvalue>::Map(omegaVector.data(), omegaVector.size()); Eigen::VectorX<fvalue> omega = Eigen::VectorX<fvalue>::Map(omegaVector.data(), omegaVector.size());
Eigen::VectorX<fvalue> drt = calcDrt<fvalue>(impedanceSpectra, omega, fm, fp); Eigen::VectorX<fvalue> drt = calcDrt<fvalue>(impedanceSpectra, omega, fm, fp, rSeries);
std::vector<fvalue> stdvector(drt.data(), drt.data()+drt.size()); std::vector<fvalue> stdvector(drt.data(), drt.data()+drt.size());
return stdvector; return stdvector;
} }
std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, FitMetics& fm, const FitParameters& fp) std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, FitMetics& fm, const FitParameters& fp, fvalue* rSeries)
{ {
Eigen::VectorX<fvalue> omega; Eigen::VectorX<fvalue> omega;
Eigen::VectorX<std::complex<fvalue>> impedanceSpectra = eistoeigen(data, &omega); Eigen::VectorX<std::complex<fvalue>> impedanceSpectra = eistoeigen(data, &omega);
Eigen::VectorX<fvalue> drt = calcDrt<fvalue>(impedanceSpectra, omega, fm, fp); Eigen::VectorX<fvalue> drt = calcDrt<fvalue>(impedanceSpectra, omega, fm, fp, rSeries);
std::cout<<omega.size()<<','<<impedanceSpectra.size()<<','<<drt.size()<<std::endl;
std::vector<fvalue> stdvector(drt.data(), drt.data()+drt.size()); std::vector<fvalue> stdvector(drt.data(), drt.data()+drt.size());
return stdvector; return stdvector;
} }
......
...@@ -37,11 +37,13 @@ Api for use with Eigen applications ...@@ -37,11 +37,13 @@ Api for use with Eigen applications
* @param omegaTensor vector with the omega values that the impedances where mesured at * @param omegaTensor vector with the omega values that the impedances where mesured at
* @param fm a fit metrics struct where this function returns information on the fit aquired * @param fm a fit metrics struct where this function returns information on the fit aquired
* @param fp a struct with fit parameters * @param fp a struct with fit parameters
* @param rSeries an optional paramter where the seires resistance is stored
* @return a vector with the drt values * @return a vector with the drt values
*/ */
template<typename fv> template<typename fv>
Eigen::VectorX<fv> calcDrt(Eigen::VectorX<std::complex<fv>>& impedanceSpectra, Eigen::VectorX<fv>& omegaTensor, FitMetics& fm, const FitParameters& fp); Eigen::VectorX<fv> calcDrt(Eigen::VectorX<std::complex<fv>>& impedanceSpectra, Eigen::VectorX<fv>& omegaTensor,
FitMetics& fm, const FitParameters& fp, fv* rSeries = nullptr);
/** /**
.... ....
......
...@@ -36,9 +36,11 @@ Api for use with eisgenerator applications ...@@ -36,9 +36,11 @@ Api for use with eisgenerator applications
* @param omegaVector vector with the omega values that the impedances where mesured at * @param omegaVector vector with the omega values that the impedances where mesured at
* @param fm a fit metrics struct where this function returns information on the fit aquired * @param fm a fit metrics struct where this function returns information on the fit aquired
* @param fp a struct with fit parameters * @param fp a struct with fit parameters
* @param rSeries an optional paramter where the seires resistance is stored
* @return a vector with the drt values * @return a vector with the drt values
*/ */
std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, const std::vector<fvalue>& omegaVector, FitMetics& fm, const FitParameters& fp); std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, const std::vector<fvalue>& omegaVector, FitMetics& fm,
const FitParameters& fp, fvalue* rSeries = nullptr);
/** /**
* @brief calculate a drt on eisgenerator types * @brief calculate a drt on eisgenerator types
...@@ -46,9 +48,10 @@ std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, const std:: ...@@ -46,9 +48,10 @@ std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, const std::
* @param data a vector of eisgenerator datapoints with the values to your expirament, embedded omega values are used * @param data a vector of eisgenerator datapoints with the values to your expirament, embedded omega values are used
* @param fm a fit metrics struct where this function returns information on the fit aquired * @param fm a fit metrics struct where this function returns information on the fit aquired
* @param fp a struct with fit parameters * @param fp a struct with fit parameters
* @param rSeries an optional paramter where the seires resistance is stored
* @return a vector with the drt values * @return a vector with the drt values
*/ */
std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, FitMetics& fm, const FitParameters& fp); std::vector<fvalue> calcDrt(const std::vector<eis::DataPoint>& data, FitMetics& fm, const FitParameters& fp, fvalue* rSeries = nullptr);
/** /**
.... ....
......
...@@ -79,12 +79,16 @@ int main(int argc, char** argv) ...@@ -79,12 +79,16 @@ int main(int argc, char** argv)
FitMetics fm; FitMetics fm;
// calculate the drt for this spectrum // calculate the drt for this spectrum
std::vector<fvalue> x = calcDrt(data, fm, FitParameters(1000)); fvalue rSeries;
std::vector<fvalue> x = calcDrt(data, fm, FitParameters(1000), &rSeries);
assert(x.size() == data.size());
// print some info on the drt // print some info on the drt
std::cout<<"Iterations: "<<fm.iterations<<'\n'; std::cout<<"Iterations: "<<fm.iterations<<'\n';
std::cout<<"fx "<<fm.fx<<"\ndrt: "; std::cout<<"fx "<<fm.fx<<"\ndrt: ";
printFvalueVector(x); printFvalueVector(x);
std::cout<<"r series: "<<rSeries<<'\n';
return 0; return 0;
} }
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment