basicmath.cpp

//SPDX-License-Identifier:         LGPL-3.0-or-later
//
// eisgenerator - a shared library and application to generate EIS spectra
// Copyright (C) 2022-2024 Carl Philipp Klemm <carl@uvos.xyz>
//
// This file is part of eisgenerator.
//
// eisgenerator is free software: you can redistribute it and/or modify
// it under the terms of the GNU Lesser General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// eisgenerator is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU Lesser General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public License
// along with eisgenerator.  If not, see <http://www.gnu.org/licenses/>.
//

#include "basicmath.h"
#include <algorithm>
#include <random>
#include <limits>
#include <stdexcept>

#include "log.h"
#include "linearregession.h"

static size_t gradIndex(size_t dataSize, size_t inputIndex)
{
	if(inputIndex == 0)
		inputIndex = 1;
	else if(inputIndex > dataSize-2)
		inputIndex = dataSize-2;
	return inputIndex;
}

std::complex<fvalue> eis::absGrad(const std::vector<eis::DataPoint>& data, size_t index)
{
	if(data.size() < 3)
		return std::complex<fvalue>(1,1);

	index = gradIndex(data.size(), index);

	return std::complex<fvalue>(std::abs((data[index+1].im.real()-data[index-1].im.real())/(data[index+1].omega-data[index-1].omega)),
								std::abs((data[index+1].im.imag()-data[index-1].im.imag())/(data[index+1].omega-data[index-1].omega)));
}

void eis::mulitplyAdd(std::vector<eis::DataPoint>& data, fvalue mult, fvalue add)
{
	for(eis::DataPoint& point : data)
		point.im = point.im*mult + add;
}

fvalue eis::grad(const std::vector<fvalue>& data, const std::vector<fvalue>& omega, size_t index)
{
	assert(data.size() == omega.size());
	if(data.size() < 3)
		return 0;

	index = gradIndex(data.size(), index);

	return (data[index+1]-data[index-1])/(omega[index+1]-omega[index-1]);
}

std::complex<fvalue> eis::grad(const std::vector<eis::DataPoint>& data, size_t index)
{
	if(data.size() < 3)