DSCAnalysisTool/src/data/pointcalculate.cpp

#include <QDebug>

#include "pointcalculate.h"
#include "logger.h"

QVector<Global::ExperimentData>PointCalculate:: _dataVtr;
QPointF PointCalculate::_peakPoint;
QPointF PointCalculate::_leftSelectedPoint,PointCalculate::_rightSelectedPoint;

void PointCalculate::setExperimentData(const QVector<Global::ExperimentData> &dataVtr)
{
    _dataVtr = dataVtr;

    if(_dataVtr.empty())    {
        return;
    }

    Global::ExperimentData startPoint = _dataVtr.at(0);
    Global::ExperimentData endPoint = _dataVtr.at(_dataVtr.size() - 1);

    _leftSelectedPoint = QPointF(startPoint.sampleTemp,startPoint.dsc);
    _rightSelectedPoint = QPointF(endPoint.sampleTemp,endPoint.dsc);
}

std::vector<float> PointCalculate::movingAverage(const std::vector<float>& data, int windowSize) {
    std::vector<float> smoothedData;
    for (size_t i = 0; i < data.size(); ++i) {
        if (i + windowSize <= data.size()) {
            float sum = std::accumulate(data.begin() + i, data.begin() + i + windowSize, 0.0);
            smoothedData.push_back(sum / windowSize);
        }
    }
    return smoothedData;
}

QPointF PointCalculate::getPeakPoint(){
    int n = _dataVtr.size();

    if (n < 3) {
        return QPointF(); // 至少需要三个点才能找到波峰
    }

    QPointF uniquePeak;
    float maxDiff = -std::numeric_limits<float>::infinity();

    for (int i = 0; i < n; ++i) {
        const float currentX = _dataVtr.at(i).sampleTemp;
        const float currentY = _dataVtr.at(i).dsc;

        if (currentX < _leftSelectedPoint.x()) {
            continue;
        }
        if (currentX > _rightSelectedPoint.x()) {
            break;
        }

        // 计算当前点与左选择点 y 值的差值
        float diffLeft = std::abs(currentY - _leftSelectedPoint.y());
        // 计算当前点与右选择点 y 值的差值
        float diffRight = std::abs(currentY - _rightSelectedPoint.y());
        // 取两个差值中的较大值
        float currentDiff = std::max(diffLeft, diffRight);

        if (currentDiff > maxDiff) {
            maxDiff = currentDiff;
            uniquePeak = QPointF(currentX, currentY);
        }
    }

    _peakPoint = uniquePeak;

    logde<<"peakPoint:"<<_peakPoint.x()<<","<<_peakPoint.y();

    return uniquePeak;
}

QPair<QPointF, QPointF> PointCalculate::calculateMaxDiffPointDetail(
        const PointCalculate::MaxDiffPointDetailType type)
{
#if 1
    float maxDiff = std::numeric_limits<float>::min();
    QPointF currentPoint,lastPoint;

    for (int i = 0; i < _dataVtr.size() - 1; ++i) {
        const float currentX = _dataVtr.at(i).sampleTemp;
        const float currentY = _dataVtr.at(i).dsc;

        if(type == MaxDiffPointDetailType::Left){
            if(currentX <= _leftSelectedPoint.x()){
                continue;
            }
            if(currentX >= _peakPoint.x()){
                break;
            }
        }else{
            if(currentX <= _peakPoint.x()){
                continue;
            }
            if(currentX >= _rightSelectedPoint.x()){
                break;
            }
        }
        //
        const float lastX = _dataVtr.at(i + 1).sampleTemp;
        const float lastY = _dataVtr.at(i + 1).dsc;
        float diff = std::abs(currentY - lastY);
        if(diff > maxDiff){
            maxDiff = diff;

            currentPoint.setX(currentX);
            currentPoint.setY(currentY);
            lastPoint.setX(lastX);
            lastPoint.setY(lastY);
        }
    }
#endif

    return qMakePair(currentPoint,lastPoint);
}


QPair<QPointF, QPointF> PointCalculate::calculateMaxDiffPointLeft()
{
    return calculateMaxDiffPointDetail(MaxDiffPointDetailType::Left);
}

QPair<QPointF, QPointF> PointCalculate::calculateMaxDiffPointRight()
{
    return calculateMaxDiffPointDetail(MaxDiffPointDetailType::Right);
}

#if 0
QPointF PointCalculate::findClosestY(float targetX)
{
    float minDiff = std::numeric_limits<float>::max();
    QPointF resultPointF;

    for(FileManager::ExperimentData &ed:_dataVtr){
        float diff = std::abs(ed.sampleTemp - targetX);
        if (diff < minDiff) {
            minDiff = diff;
            resultPointF = QPointF(ed.sampleTemp,ed.dsc);
        }
    }

    return resultPointF;
}
#endif

#if 0
QPointF PointCalculate::findClosestPointByX(float x) {
    int left = 0;
    int right = _dataVtr.size() - 1;

    QPointF targetPoint;
    targetPoint.setX(_dataVtr.value(0).sampleTemp);
    targetPoint.setY(_dataVtr.value(0).dsc);

    while (left <= right) {
        int mid = left + (right - left) / 2;
        FileManager::ExperimentData& ed = _dataVtr[mid];

        if (std::abs(ed.sampleTemp - x) < std::abs(targetPoint.x() - x)) {
            targetPoint.setX(ed.sampleTemp);
            targetPoint.setY(ed.dsc);
        }

        if(ed.sampleTemp < x){
            left = mid + 1;
        } else {
            right = mid - 1;
        }
    }

    return targetPoint;
}
#endif

void PointCalculate::setRegionPointX(const float left, const float right)
{
    logde<<"dataVtr size:"<<_dataVtr.size();
    logde<<"select point param left,right:"<<left<<","<<right;

    _leftSelectedPoint = getClosestPointByX(left);
    _rightSelectedPoint = getClosestPointByX(right);

    _peakPoint = getPeakPoint();

    logde<<"peak point:"<<_peakPoint.x()<<","<<_peakPoint.y();

    //根据峰谷重置选择点。
    updateStartEndPoint();

    logde<<"select point left:"<<_leftSelectedPoint.x()<<","<<_leftSelectedPoint.y();
    logde<<"select point right:"<<_rightSelectedPoint.x()<<","<<_rightSelectedPoint.y();
}

QVector<QPointF> PointCalculate::getPeakPointGroup()
{
    QVector<QPointF> pointVtr;
    for(Global::ExperimentData& ed:_dataVtr)   {
        if(ed.sampleTemp >= _leftSelectedPoint.x() &&
                ed.sampleTemp <= _rightSelectedPoint.x()){
            pointVtr.push_back(QPointF(ed.sampleTemp,ed.dsc));
        }
    }

    return pointVtr;
}

float PointCalculate::calculateArea() {
    //getPoint group
    QVector<QPointF> points = getPeakPointGroup();

    //calculate Area
    float integral = 0.0;
    size_t n = points.size();

    if (n < 2) {
        return integral; // 至少需要两个点才能计算积分
    }

    //find a line.
    float k = (_leftSelectedPoint.y() - _rightSelectedPoint.y()) /
            (_leftSelectedPoint.x() - _rightSelectedPoint.x());
    float b = _leftSelectedPoint.y() - k * _leftSelectedPoint.x();

    for (size_t i = 0; i < n - 1; ++i) {
#if 1
        float x1 = points[i].x();
        float y1 = points[i].y();

        float x2 = points[i + 1].x();
        float y2 = points[i + 1].y();
#endif
        float yLine1 = k * x1 + b;
        float yLine2  = k * x2 + b;

        float diff1 = y1 - yLine1;
        float diff2 = y2 - yLine2;

        float dx = x2 - x1;

        float cellArea = (diff1 + diff2) * dx /2.0;

        integral += std::abs(cellArea);
    }


    return integral;
}

QPair<QPointF, QPointF> PointCalculate::calculateStartAndEndPoint()
{
    QPair<QPointF,QPointF> leftMaxDiffPointPair =  PointCalculate::calculateMaxDiffPointLeft();
    QPair<QPointF,QPointF> rightMaxDiffPointPair =  PointCalculate::calculateMaxDiffPointRight();
#if 0
    logde<<"b1:"<<leftMaxDiffPointPair.first.x()<<","<<leftMaxDiffPointPair.first.y();
    logde<<"b2:"<<leftMaxDiffPointPair.second.x()<<","<<leftMaxDiffPointPair.second.y();

    logde<<"b3:"<<rightMaxDiffPointPair.first.x()<<","<<rightMaxDiffPointPair.first.y();
    logde<<"b4:"<<rightMaxDiffPointPair.second.x()<<","<<rightMaxDiffPointPair.second.y();
#endif
    QPointF startPoint = calculateIntersection(_leftSelectedPoint,_rightSelectedPoint,
                                               leftMaxDiffPointPair.first,
                                               leftMaxDiffPointPair.second);

    QPointF endPoint = calculateIntersection(_leftSelectedPoint,_rightSelectedPoint,
                                             rightMaxDiffPointPair.first,
                                             rightMaxDiffPointPair.second);

    return qMakePair(startPoint,endPoint);
}
void PointCalculate::updateStartEndPoint()
{
    //需要在a1和a2之间查询是否有高于a1和a2之间的点，若存在，则重新给a1、a2赋值

    for(Global::ExperimentData& ed:_dataVtr){
        if(ed.sampleTemp > _leftSelectedPoint.x() && ed.sampleTemp < _peakPoint.x()){
            if(ed.dsc > _leftSelectedPoint.y()){
                _leftSelectedPoint.setX(ed.sampleTemp);
                _leftSelectedPoint.setY(ed.dsc);
            }
        }

        if(ed.sampleTemp < _rightSelectedPoint.x() && ed.sampleTemp > _peakPoint.x()){
            if(ed.dsc > _rightSelectedPoint.y()){
                _rightSelectedPoint.setX(ed.sampleTemp);
                _rightSelectedPoint.setY(ed.dsc);
            }
        }
    }
}
QString PointCalculate::textFormatPeakPoint(const float enthalpyValue,
                                            const float peakValue,
                                            const float startPoint,
                                            const float endPoint)
{
    return QString("峰的综合信息：\n"
                   "焓值：%1 J/g \n"
                   "峰值：%2℃ \n"
                   "起始点：%3℃ \n"
                   "终止点：%4℃"
                   ).arg(QString::number(enthalpyValue, 'f', 3))
            .arg(QString::number(peakValue, 'f', 3))
            .arg(QString::number(startPoint, 'f', 3))
            .arg(QString::number(endPoint, 'f', 3));
}

// 计算两条直线的交点
QPointF PointCalculate::calculateIntersection(const QPointF p1,const QPointF p2,
                                              const QPointF p3, const QPointF p4){
    // 直线的一般式: A1x + B1y + C1 = 0 和 A2x + B2y + C2 = 0
    float A1 = p2.y() - p1.y();
    float B1 = p1.x() - p2.x();
    float C1 = A1 * p1.x() + B1 * p1.y();

    float A2 = p4.y() - p3.y();
    float B2 = p3.x() - p4.x();
    float C2 = A2 * p3.x() + B2 * p3.y();

    float determinant = A1 * B2 - A2 * B1;

    if (determinant == 0) {
        // 两条直线平行或重合，无交点或无限交点
        return {0, 0};
    } else {
        float x = (B2 * C1 - B1 * C2) / determinant;
        float y = (A1 * C2 - A2 * C1) / determinant;
        return {x, y};
    }
}

QPair<float,float> PointCalculate::getCurveInflectionPointTangent(const float x1,const float x2)
{
    std::vector<float> dataVtr;
    for(Global::ExperimentData& ed:_dataVtr){
        if(x1 < ed.sampleTemp && ed.sampleTemp < x2){
            dataVtr.push_back(ed.sampleTemp);
        }
    }
    //
    std::vector<float> processedVtr = movingAverage(dataVtr,5);

    // find max slope.
    float maxSlope = 0.0;
    float targetX = 0.0;
    for (size_t i = 1; i < processedVtr.size(); ++i) {
        float slope = processedVtr[i] - processedVtr[i - 1];
        if(std::abs(slope) > maxSlope){
            maxSlope = slope;
            targetX = processedVtr[i];
        }
    }

    // calculate line slope and axis Y space number b.
    QPointF point = getClosestPointByX(targetX);

    // 使用点斜式方程求y轴截距
    // y = m * x + b
    //      double b = y - m * x;
    float b = point.y() - maxSlope * point.x();

    return qMakePair<float,float>(maxSlope,b);
}
QPointF PointCalculate::getIntersectionBySlope(const LineStruct& line1, const LineStruct& line2) {
    float x = (line2.intercept - line1.intercept) / (line1.slope - line2.slope);
    float y = line1.slope * x + line1.intercept;
    //    return {x, y};
    return QPointF(x,y);
}
QPointF PointCalculate::getClosestPointByX(const float targetX)
{
#if 1
    QPointF resultPointF;
    float minDiff = std::numeric_limits<float>::max();

    for(Global::ExperimentData &ed:_dataVtr){
        float diff = std::abs(ed.sampleTemp - targetX);
        if (diff < minDiff) {
            minDiff = diff;
            resultPointF = QPointF(ed.sampleTemp,ed.dsc);
        }
    }

    return resultPointF;
#endif

#if 0
    int left = 0;
    int right = _dataVtr.size() - 1;

    while (left < right) {
        int mid = left + (right - left)/2;
        const FileManager::ExperimentData& ed = _dataVtr.at(mid);
        if(ed.sampleTemp > targetX){
            right = mid;
        }else{
            left = mid + 1;
        }
    }

    //
    QPointF resultPointF(_dataVtr.at(left).sampleTemp,
                         _dataVtr.at(left).dsc);

    if(left > 0 &&
            std::abs(resultPointF.x() - targetX) >
            std::abs(_dataVtr.at(left - 1).sampleTemp - targetX)){
        //
        resultPointF = QPointF(_dataVtr.at(left - 1).sampleTemp,
                               _dataVtr.at(left - 1).dsc);
    }

    return resultPointF;
#endif
}

QString PointCalculate::textFormatNumbericalLabel(const QPointF point)
{
    return QString("数值：\n"
                   "%1℃,%2"
                   ).arg(QString::number(point.x(), 'f', 3))
            .arg(QString::number(point.y(), 'f', 3));
}

QPair<QPointF, QPointF> PointCalculate::getStartAndEndPoint()
{
    if(_dataVtr.empty()){
        return qMakePair(QPointF(), QPointF());
    }

    Global::ExperimentData startPoint  = _dataVtr.at(0);
    Global::ExperimentData endPoint    = _dataVtr.at(_dataVtr.size() - 1);

    return qMakePair<QPointF,QPointF>(
                QPointF(startPoint.sampleTemp,startPoint.dsc),
                QPointF(endPoint.sampleTemp,endPoint.dsc));
}

QString PointCalculate::textFormatStartPoint(const QPointF point)
{
    return QString("外推起始点：\n"
                   "%1℃"
                   ).arg(QString::number(point.x(), 'f', 3));
}

QString PointCalculate::textFormatEndPoint(const QPointF point)
{
    return QString("外推终止点：\n"
                   "%1℃"
                   ).arg(QString::number(point.x(), 'f', 3));
}


QPointF PointCalculate::getClosestPointByY(const float left,const float right,const float valueY)
{
    float minValue = std::numeric_limits<float>::infinity(); // 初始化为正无穷
    QPointF closestPoint;

    for (const Global::ExperimentData& ed : _dataVtr) {
        if (left < ed.sampleTemp && ed.sampleTemp < right) {
            float diff = std::abs(ed.dsc - valueY);
            if (diff < minValue) {
                minValue = diff;
                closestPoint = QPointF(ed.sampleTemp, ed.dsc);
            }
        }
    }

    return closestPoint;
}

QString PointCalculate::textFormatGlassTranstion(const float t1,const float tg,const float t2)
{
    return QString("T1：%1℃\n"
                   "Tg：%2℃\n"
                   "T2：%3℃"
                   ).arg(QString::number(t1, 'f', 3))
            .arg(QString::number(tg, 'f', 3))
            .arg(QString::number(t2, 'f', 3));
}


QPair<float, float> PointCalculate::getMaxMinValue()
{
    if (_dataVtr.isEmpty()) {
        return QPair<float, float>(0.0f, 0.0f);
    }

    float maxDsc = std::numeric_limits<float>::lowest(); // 初始化为最小浮点数
    float minDsc = std::numeric_limits<float>::max();    // 初始化为最大浮点数

    for(const Global::ExperimentData &ed : _dataVtr) {
        if (ed.dsc > maxDsc) {
            maxDsc = ed.dsc;
        }
        if (ed.dsc < minDsc) {
            minDsc = ed.dsc;
        }
    }

    return QPair<float, float>(minDsc, maxDsc);
}

QVector<Global::ExperimentData> PointCalculate::getDataInXRange(const float x1, const float x2)
{
    QVector<Global::ExperimentData> targetVtr;
    for(const Global::ExperimentData &ed : _dataVtr) {
        if(x1 < ed.sampleTemp && ed.sampleTemp < x2){
            targetVtr.push_back(ed);
        }else if (ed.sampleTemp > x2){
            break;
        }

    }

    return targetVtr;
}