DSCAnalysisTool/src/data/peakpoint.cpp

#include <QDebug>

#include "peakpoint.h"
#include "logger.h"

QVector<FileManager::ExperimentData>PeakPoint:: _dataVtr;
QPointF PeakPoint::_peakPoint;
QPointF PeakPoint::_leftSelectedPoint,PeakPoint::_rightSelectedPoint;

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

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

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

QPointF PeakPoint::findPeakPoint(){
    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> PeakPoint::calculateMaxDiffPointDetail(
        const PeakPoint::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> PeakPoint::calculateMaxDiffPointLeft()
{
    return calculateMaxDiffPointDetail(MaxDiffPointDetailType::Left);
}

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

#if 0
QPointF PeakPoint::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 PeakPoint::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 PeakPoint::setRegionPointX(const float left, const float right)
{
    logde<<"select point param left,right:"<<left<<","<<right;

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

    _peakPoint = findPeakPoint();
    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> PeakPoint::getPeakPointGroup()
{
    QVector<QPointF> pointVtr;
    for(FileManager::ExperimentData& ed:_dataVtr)   {
        if(ed.sampleTemp >= _leftSelectedPoint.x() && ed.sampleTemp <= _rightSelectedPoint.x()){
            pointVtr.push_back(QPointF(ed.sampleTemp,ed.dsc));
        }
    }

    return pointVtr;
}

double PeakPoint::calculateArea() {
    //getPoint group
    QVector<QPointF> points = getPeakPointGroup();
#if 0
    for(QPointF point:points){
        qDebug()<<"x,y:"<<QString::number(point.x(),'f',3)<<","
               <<QString::number(point.y(),'f',3);
    }
#endif
#if 0
    qDebug()<<QString::number(points.at(0).x(),'f',3)<<","<<
              QString::number(points.at(points.size() - 1).x(),'f',3);
#endif

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

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

#if 0
    for (size_t i = 1; i < n; ++i) {
        float dx = points[i].x() - points[i - 1].x();
        float avg_y = (points[i].y() + points[i - 1].y()) / 2.0;

#if 0
        float step =  dx *  std::abs(avg_y);
        qDebug()<<points[i].x()<<","<<points[i].y()<<"integral:"<<integral<<","<<step;
#endif

        integral += dx * std::abs(avg_y);
    }
#endif
#if 0
    double k = (y2 - y1) / (x2 - x1);
    double b = y1 - k * x1;
#endif
    //find a line.
    double k = (_leftSelectedPoint.y() - _rightSelectedPoint.y()) /
            (_leftSelectedPoint.x() - _rightSelectedPoint.x());
    double b = _leftSelectedPoint.y() - k * _leftSelectedPoint.x();

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

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

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

        double dx = x2 - x1;

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

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


    return integral;
}

QPair<QPointF, QPointF> PeakPoint::calculateStartAndEndPoint()
{
    QPair<QPointF,QPointF> leftMaxDiffPointPair =  PeakPoint::calculateMaxDiffPointLeft();
    QPair<QPointF,QPointF> rightMaxDiffPointPair =  PeakPoint::calculateMaxDiffPointRight();

    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();

    QPointF startPoint = calculateIntersection(_leftSelectedPoint,_rightSelectedPoint,
                                               leftMaxDiffPointPair.first,
                                               leftMaxDiffPointPair.second);

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

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

    for(FileManager::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 PeakPoint::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 PeakPoint::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};
    }
}

QPointF PeakPoint::findClosestPointByX(const 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;
}

QString PeakPoint::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> PeakPoint::getStartAndEndPoint()
{
    FileManager::ExperimentData startPoint = _dataVtr.at(0);
    FileManager::ExperimentData endPoint = _dataVtr.at(_dataVtr.size() - 1);

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