用hist函数计算直方图的MATLAB算法是什么?

ldioqlga  于 2022-11-24  发布在  Matlab
关注(0)|答案(1)|浏览(143)

我正在将一些旧的MATLAB代码翻译成C++。我注意到,我的自定义函数计算直方图,它应该等效于MATLAB [counts,centers]= hist(___),但却给出了不同的结果。我在我的实现中找不到错误。所以我使用MATLAB编码器从MATLAB代码生成C函数,并将其与我的C代码进行比较。密码:

function [counts, centers] = my_hist(values, bins)
    [counts, centers] = hist(values, bins);
    disp(centers);
    disp(counts);
end

和一个调用它的脚本,这样MATLAB就可以定义输入:

values = rand(1,1000);
bins = linspace(0.05, 0.95, 10);

[counts, centers] = my_hist(values, bins);

基于以上,编码器生成函数:

//
// File: my_hist.cpp
//
// MATLAB Coder version            : 5.3
// C/C++ source code generated on  : 17-Nov-2022 15:46:17
//

// Include Files
#include "my_hist.h"
#include "rt_nonfinite.h"
#include <algorithm>
#include <cmath>
#include <cstring>
#include <math.h>

// Function Definitions
//
// MY_HIST Summary of this function goes here
//    Detailed explanation goes here
//
// Arguments    : const double values[1000]
//                const double bins[10]
//                double counts[10]
//                double centers[10]
// Return Type  : void
//
void my_hist(const double values[1000], const double bins[10],
             double counts[10], double centers[10])
{
  double edges[11];
  double nn[11];
  double absx;
  int k;
  int low_i;
  std::copy(&bins[0], &bins[10], &centers[0]);
  for (k = 0; k < 9; k++) {
    absx = bins[k];
    edges[k + 1] = absx + (bins[k + 1] - absx) / 2.0;
  }
  edges[0] = rtMinusInf;
  edges[10] = rtInf;
  for (k = 0; k < 9; k++) {
    double absx_tmp;
    absx_tmp = edges[k + 1];
    absx = std::abs(absx_tmp);
    if ((!std::isinf(absx)) && (!std::isnan(absx))) {
      if (absx <= 2.2250738585072014E-308) {
        absx = 4.94065645841247E-324;
      } else {
        frexp(absx, &low_i);
        absx = std::ldexp(1.0, low_i - 53);
      }
    } else {
      absx = rtNaN;
    }
    edges[k + 1] = absx_tmp + absx;
  }
  std::memset(&nn[0], 0, 11U * sizeof(double));
  low_i = 1;
  int exitg1;
  do {
    exitg1 = 0;
    if (low_i + 1 < 12) {
      if (!(edges[low_i] >= edges[low_i - 1])) {
        for (low_i = 0; low_i < 11; low_i++) {
          nn[low_i] = rtNaN;
        }
        exitg1 = 1;
      } else {
        low_i++;
      }
    } else {
      for (k = 0; k < 1000; k++) {
        low_i = 0;
        absx = values[k];
        if (!std::isnan(absx)) {
          if ((absx >= edges[0]) && (absx < edges[10])) {
            int high_i;
            int low_ip1;
            low_i = 1;
            low_ip1 = 2;
            high_i = 11;
            while (high_i > low_ip1) {
              int mid_i;
              mid_i = (low_i + high_i) >> 1;
              if (values[k] >= edges[mid_i - 1]) {
                low_i = mid_i;
                low_ip1 = mid_i + 1;
              } else {
                high_i = mid_i;
              }
            }
          }
          if (values[k] == edges[10]) {
            low_i = 11;
          }
        }
        if (low_i > 0) {
          nn[low_i - 1]++;
        }
      }
      exitg1 = 1;
    }
  } while (exitg1 == 0);
  std::copy(&nn[0], &nn[10], &counts[0]);
  counts[9] += nn[10];
}

//
// File trailer for my_hist.cpp
//
// [EOF]
//

我不明白在这段代码中发生了什么以及为什么这样做:

for (k = 0; k < 9; k++) {
    double absx_tmp;
    absx_tmp = edges[k + 1];
    absx = std::abs(absx_tmp);
    if ((!std::isinf(absx)) && (!std::isnan(absx))) {
      if (absx <= 2.2250738585072014E-308) {
        absx = 4.94065645841247E-324;
      } else {
        frexp(absx, &low_i);
        absx = std::ldexp(1.0, low_i - 53);
      }
    } else {
      absx = rtNaN;
    }
    edges[k + 1] = absx_tmp + absx;
  }

该功能移动垃圾箱的边缘,但如何和为什么?我将感谢帮助和解释!

lymnna71

lymnna711#

这段代码将eps添加到除第一个和最后一个bin边缘之外的每个bin边缘。
很难知道hist为什么这样做,他们一定是在处理他们发现的一些边缘情况(可能与浮点舍入错误有关),并认为这是最好或最简单的解决方案。

相关问题