2023年第29周个人总结

发布时间 2023-07-15 13:25:51作者: akaS

文档编辑机器+移动窗口变化率计算

  1. 投标文档书写

    本周在公司主要对投标文档的修改:

    • 首先表格的脚注很重要,当多个表格在一起时,表格的脚注可以起到非常关键的作用,能够让读者快速弄懂该表格的含义。
    • 然后关于文档的错别字,自己写的文档自己是真的不可能能够读懂,所以可以采用通读一遍的方式进行纠错。
    • 关于表格的断表(即因为分页后导致的表格不连续),需要在表头的地方加一个表格连续,表格属性-行-最后一个选项。
  2. python代码移植c++

    本周python程序移植主要做了关于我们对方位数据分窗口采用最小二乘法平滑进行变化率的求解,因为我们最开始的idea无法实现对带噪声的误差数据进行求方位变化率,会导致失真。所以我们采用窗口法,先定一个窗口大小len_,然后确定一个窗口移动距离cnt,再根据数据总长度结合前两个参数,可以确定对该段数据需要做多少次的窗口计算。

#最小二乘拟合函数
def err_func(p,x,y):
    return fit_func(p,x)-y

def my_leastsq(err_func, x,y):
    p0=[0.01,-0.02,0.03,-0.02]
    params, success = optimize.leastsq(err_func, p0, args=(x,y))
    y_smooth = fit_func(params,x)
    return y_smooth

def curve_smooth(y):
    x = np.arange(len(y)).astype(np.float32)
    y_smooth = my_leastsq(err_func,x,y)
    return y_smooth

#滑动窗口法计算变化率
def resize_list(y,len_y):
    interp_fn = interpld(np.arrange(len(y)).astype(np.float32),y,kind='cubic')
    y_ = interp_fn(np.linspace(0,len(y)-1,len_y_))
        return y_

def area_rate(ed,ei_center,si_center,radius):
    d = np.abs(ed[ei_center-radius:ei_center+radius-1]-ed[si_center-radius:si_center+radius-1])
    return np.median(d)/(ei_center-si_center)*60
    
def brg_rate_curve(e_n):
    len_ = 60
    fwd_steps = len_//2
    step_cnt = (len(e_n)-len)//(fwd_steps)-1
    rate_set = np.zeros(step_cnt+1,np.float32)
    for i in range(step_cnt)
        si = i * fwd_steps
        ei = si + len_
        edata_n = e_n[si:ei]
        edata_s_ = curve_smooth(edata_n)
        sam_ei = len_-len_//4
        sam_si = len_//4
        rate_set[i+1] = area_rate(edata_s_,sam_ei,sam_si,3)#分段进行窗口计算
        pass
    rate_set = np.abs(rate_set)
    rate_set[0] = (rate_set[1] - 0.50092*rate_set[2])/(1-0.50092)
    rate_set_sz = resize_list(rate_set,len(e_n))#返回插值结果,将数据结果进行平滑补点
    e_s_ = curve_smooth(e_n)
    brg_c_ = e_s_
    brg_rate_c = rate_set_sz
    
    return brg_c_ , brg_rate_c

​ 通过将上述python代码移植到C++环境中进行实现,由于实际应用场景的编译器首先,不可以使用各种外接库,故需要从数学原理上对上述逻辑进行实现。

#include <iostream>
#include<fstream>
#include <vector>
#include <cmath>
#include <vector>
#include <utility>
#include <sstream>
#include <cmath>
#include <algorithm>

std::vector<double> leastSquaresSmoothing_1in(const std::vector<double>& coordinates, int degree)
{
    int n = coordinates.size();
    std::vector<double> x(n);
    std::vector<double> y(n);
    for (int i = 0; i < n; ++i) {
        x[i] = i;
        y[i] = coordinates[i];
    }
    // 构建矩阵 A 和向量 b
    std::vector<std::vector<double>> A(degree + 1, std::vector<double>(degree + 1, 0.0));
    std::vector<double> b(degree + 1, 0.0);
    for (int i = 0; i <= degree; ++i) {
        for (int j = 0; j <= degree; ++j) {
            for (int k = 0; k < n; ++k) {
                A[i][j] += pow(x[k], i + j);
            }
        }
    }
    for (int i = 0; i <= degree; ++i) {
        for (int k = 0; k < n; ++k) {
            b[i] += pow(x[k], i) * y[k];
        }
    }
    // 解线性方程组
    for (int i = 0; i <= degree; ++i) {
        for (int j = i + 1; j <= degree; ++j) {
            double ratio = A[j][i] / A[i][i];
            for (int k = i; k <= degree; ++k) {
                A[j][k] -= ratio * A[i][k];
            }
            b[j] -= ratio * b[i];
        }
    }

    for (int i = degree; i >= 0; --i) {
        b[i] /= A[i][i];
        for (int j = i - 1; j >= 0; --j) {
            b[j] -= A[j][i] * b[i];
        }
    }

    // 计算平滑后的坐标数据
    std::vector<double> smoothedCoordinates(n);
    for (int i = 0; i < n; ++i) {
        double smoothedY = 0.0;
        for (int j = 0; j <= degree; ++j) {
            smoothedY += b[j] * pow(x[i], j);
        }
        smoothedCoordinates[i] = smoothedY;
    }

    return smoothedCoordinates;
}

// Function to perform cubic interpolation
//std::vector<double> interpolate(const std::vector<double>& y, int len_y) {
//    std::vector<double> interpolated(len_y);
//    int n = y.size();
//
//    for (int i = 0; i < len_y; i++) {
//        double t = static_cast<double>(i) * (n - 1) / (len_y - 1);
//        int index = static_cast<int>(t);
//        double fractional = t - index;
//
//        if (index >= n - 1) {
//            interpolated[i] = y[n - 1];
//        }
//        else {
//            double p0 = y[index];
//            double p1 = y[index + 1];
//            double m0 = (index == 0) ? (y[1] - y[0]) : (y[index] - y[index - 1]);
//            double m1 = (index >= n - 2) ? (y[n - 1] - y[n - 2]) : (y[index + 1] - y[index]);
//
//            double a = -2 * p0 + 2 * p1 + m0 + m1;
//            double b = 3 * p0 - 3 * p1 - 2 * m0 - m1;
//            double c = m0;
//            double d = p0;
//
//            interpolated[i] = ((a * fractional + b) * fractional + c) * fractional + d;
//        }
//    }
//
//    return interpolated;
//}
std::vector<double> resize_list(std::vector<double>& y, int len_y_) {
    std::vector<double> y_(len_y_);
    int len = y.size();
    double step = (double)(len - 1) / (len_y_ - 1);
    for (int i = 0; i < len_y_; i++) {
        double idx = step * i;
        int idx1 = (int)idx;
        int idx2 = (idx1 < len - 1) ? idx1 + 1 : idx1;
        double w = idx - idx1;
        y_[i] = (1 - w) * y[idx1] + w * y[idx2];
    }
    return y_;
}

// Function to calculate the area rate
double calculateAreaRate(const std::vector<double>& ed, int ei_center, int si_center, int radius) {
    std::vector<double> d(radius * 2 - 1);
    for (int i = 0; i < radius * 2 - 1; i++) {
        d[i] = std::abs(ed[ei_center - radius + i] - ed[si_center - radius + i]);
    }
    double median = *std::max_element(d.begin(), d.end());
    return median / (ei_center - si_center) * 60;
}


// Function to calculate brg_c_ and brg_rate_c
std::pair<std::vector<double>, std::vector<double>> calculateBrgRateCurve(const std::vector<double>& e_n) {
    int len_ = 60;
    int fwd_steps = len_ / 20;
    int num_steps = (e_n.size() - len_) / fwd_steps - 1;
    std::vector<double> rate_set(num_steps + 1);
    for (int i = 0; i < num_steps; i++) {
        int si = i * fwd_steps;
        int ei = si + len_;
        std::vector<double> edata_n(e_n.begin() + si, e_n.begin() + ei);
        std::vector<double> edata_s_ = leastSquaresSmoothing_1in(edata_n,2);
        int sam_ei = len_ - len_ / 4;
        int sam_si = len_ / 4;
        rate_set[i + 1] = calculateAreaRate(edata_s_, sam_ei, sam_si, 3);
    }
    std::transform(rate_set.begin(), rate_set.end(), rate_set.begin(), [](double rate) {
        return std::abs(rate);
        });
    rate_set[0] = (rate_set[1] - 0.50092 * rate_set[2]) / (1 - 0.50092);
    std::vector<double> rate_set_sz = resize_list(rate_set, e_n.size());
    std::vector<double> e_s_ = leastSquaresSmoothing_1in(e_n,2);
    return std::make_pair(e_s_, rate_set_sz);
}


int main() {
    std::ifstream file("edata_n.txt"); // 文件见下文
    if (file.is_open())
    {
        std::vector<double> edata_n; 
        std::string line;
        while (std::getline(file, line)) {
            std::istringstream iss(line); 
            double value;
            if (iss >> value) {
                edata_n.push_back(value); 
            }
        }
        file.close();
        std::cout << "Success to open the file." << std::endl;


        std::pair<std::vector<double>, std::vector<double>> result = calculateBrgRateCurve(edata_n);
        std::vector<double> brg_c_ = result.first;
        std::vector<double> brg_rate_c = result.second;
        std::ofstream outputFile("brg_rate_c_.txt");
        std::cout << "Success to open the brg_rate_c." << std::endl;
        if (outputFile.is_open())
        {
            for (int i = 0; i < brg_rate_c.size(); ++i)
            {
                outputFile << brg_rate_c[i] << "\n";
            }
            outputFile.close();
            std::cout << "brg_rate_c data has been written to the file 'brg_rate_c.txt'." << std::endl;
        }
        else
        {
            std::cout << "Unable to open the output file." << std::endl;
        }
    }
    else
    {
        std::cout << "Failed to open the file." << std::endl;
    }

    return 0;
}



最小二乘数学原理简单介绍:

​ 最小二乘法是一种数学优化方法,用于在给定的数据集中,寻找一个模型使得该函数的预测值与数据点实际值之间的残差平方和最小。在统计学中,最小二乘法通常被用来拟合线性回归模型,但是它也可以用于非线性回归、多项式回归、傅里叶级数拟合等其它模型的拟合。

​ 最小二乘法的原理可以通过数学为公式来表示。设有一个数据集

\[(x1,y1),(x2,y2)...(xi,yi) \]

想要拟合一个函数

\[f(x)=ax**2+bx+c \]

,使得预测值与真实值之间的误差最小。可以定义误差函数:

\[e=(yi-f(xi))^2 \]

其中y是数据集的大小。

288.0722321112773
288.4733267973866
288.5553904078356
288.49691338218184
288.18874229048436
287.588055187814
287.9251878878241
287.4197593064347
288.06652695750415
288.28976700763127
287.6274003215959
288.10785911378395
287.44683125720564
288.05653436281267
287.9216680441517
287.9051569509669
286.9277872297338
287.76301062814423
286.8025457595841
287.5812292913938
286.7999670442961
287.06887490344957
287.52040888797967
287.0537525455889
287.41681690028605
286.69803181188774
286.3805704613421
286.374543477965
286.7119305660729
286.51010820826764
286.5913170276174
286.9555472250567
286.74756444131805
286.38040084782216
286.0382980496535
286.6244152873652
285.8982058221145
286.4951770259197
286.0899971946424
285.93084806074665
285.89178036405
285.82045480365525
286.4427044572006
285.7922860823515
286.33400810029866
285.77214342487287
285.3348803175391
285.8855545549703
285.87990019496215
285.84373722289996
285.71838039769915
284.97478692092426
285.68308342242716
284.98420069212204
284.893133563312
285.3686038680218
285.02151241352817
284.7020094803976
284.84542276437327
285.3710581109227
285.18408068641855
284.8854891405613
285.0970788864539
284.29719562906894
284.55851790544995
284.3268747463519
284.45345788866524
284.4868325070452
283.9842541038891
284.14453207364306
283.9116846100682
284.1572509295454
283.8777455064121
284.6730280810886
283.9555087941477
284.6300956675856
284.5660779798638
284.0430246745522
283.91912870423215
284.0902814359195
283.67856923406003
284.12520507088016
283.37112765298974
283.24874454454994
283.27181352660716
283.09052742143444
283.1532402080655
283.19740199310166
283.53787384753156
283.4601045163742
283.4763639679476
282.86655599583423
283.6573073092819
282.91448640653675
283.43116059132115
282.9686079243279
283.457339136788
282.5835999390826
282.3603284385869
283.04805498149165
283.1929975080427
283.1514836538302
282.2168069148275
282.1168192195161
282.4323796241396
282.62552483932274
282.6524756568307
282.2891564369958
282.48192368089263
282.61811407070365
282.23878890780105
282.51110856527714
282.58389578906764
281.87541060924883
282.0839614012223
282.1449668835848
282.14380481988985
281.37344178595527
281.5194013493801
281.5016724267816