主页 / user-1008405

euraad's questions

Martin Hope
euraad
Asked: 2024-01-02 17:55:23 +0800 CST
  • 6

此 C 代码执行 DBSCAN - 对带有噪声的应用程序进行基于密度的空间聚类。它是一种聚类算法，用于将无监督（无标签）数据转变为有监督（标签，例如类号）数据。而且速度也很快。DBSCAN 的缺点是消耗大量内存，但速度非常出色！

所以我所做的是创建了一个使用递归的 DBSCAN 算法。完成此操作后，我注意到如果行参数很大，则递归会导致堆栈溢出。

问题：

  1. 我想确保 C 中的递归函数必须分配在堆上，而不是堆栈上。那怎么办呢？
  2. 有什么办法可以确保这个函数 recursion_clustering 不会增加任何内存并且只使用静态内存？我尝试使用尽可能多的指针来防止内存重复。

代码：

/* Private function */
static bool recursion_clustering(const int32_t i, const size_t* row, const float* epsilon, const size_t* min_pts, const float C[], bool visited[], size_t idx[], int32_t* id, int32_t* count);

 /*
  * Create ID's of clusters
  * X[m*n]
  * idx[m] = Contains cluster ID. Noise ID is 0
  * epsilon = Raduis of the clusters
  * min_pts = Minimum points of a valid cluster
  */
void dbscan(const float X[], size_t idx[], const float epsilon, const size_t min_pts, const size_t row, const size_t column) {
    /* Create idx */
    memset(idx, 0, row * sizeof(size_t));

    /* Create pdist2 C */
    float* C = (float*)malloc(row * row * sizeof(float));
    pdist2(X, X, C, row, column, row, PDIST2_METRIC_EUCLIDEAN);

    /* Flags */
    bool* visited = (bool*)malloc(row * sizeof(bool));
    memset(visited, false, row * sizeof(bool));
    
    /* Iterate */
    int32_t i, count, id = 1;
    for (i = 0; i < row; i++) {
        if (recursion_clustering(i, &row, &epsilon, &min_pts, C, visited, idx, &id, &count)) {
            /* Set ID value */
            idx[i] = id;
            id++;
        }
    }

    /* Free */
    free(C);
    free(visited);
}

static bool recursion_clustering(const int32_t i, const size_t* row, const float* epsilon, const size_t* min_pts, const float C[], bool visited[], size_t idx[], int32_t* id, int32_t* count) {
    /* Declare variables */
    int32_t j;

    /* Check if already visited - else mark as visited */
    if (visited[i]) {
        return false;
    }
    else {
        visited[i] = true;
    }

    /* Begin first with the row of C */
    *count = 0;
    for (j = 0; j < *row; j++) {
        /* Check how many index exceeds min_pts */
        if (C[i * (*row) + j] <= *epsilon) {
            (*count)++;
        }
    }

    /* Check if count is less than min_pts */
    if (*count < *min_pts) {
        return false;
    }

    /* Iterate forward */
    for (j = i + 1; j < *row; j++) {
        if (!visited[j] && C[i * (*row) + j] <= *epsilon) {
            /* Recursion */
            if (recursion_clustering(j, row, epsilon, min_pts, C, visited, idx, id, count)){
                /* Set ID value */
                idx[j] = *id;
            }
        }
    }

    /* Iterate backward */
    for (j = i - 1; j >= 0; j--) {
        if (!visited[j]) {
            /* Recursion */
            if (recursion_clustering(j, row, epsilon, min_pts, C, visited, idx, id, count)) {
                /* Set ID value */
                idx[j] = *id;
            }
        }
    }

    /* Return true */
    return true;
}

如果您想运行此代码，您需要这个工作示例。

#define row 65
#define column 3

int main() {
    /* Define X, idx, epsilon and min_pts */
    float X[row * column] = { -0.251521, 1.045117, -1.281658,
        -1.974109, 0.278170, -1.023392,
        -0.957729, -0.977450, 0.477872,
        -0.449159, -1.016680, 0.095610,
        -1.785787, -1.403543, 0.483454,
        1.366889, -0.762590, -1.162454,
        2.129839, 0.358568, -2.118250,
        0.751071, -1.766582, 0.178434,
        -1.980847, -1.320933, -0.457778,
        -0.478030, 0.606917, -1.630624,
        3.674916, 0.088957, 0.877373,
        0.637213, 0.079176, 0.342038,
        1.142329, 0.629997, 0.311134,
        -0.878974, 0.042527, 0.736522,
        1.751637, -1.434299, -1.325140,
        1.110682, 1.091970, 1.434869,
        -0.504482, -2.504821, -1.245315,
        -0.102915, -0.203266, -0.849767,
        -0.822834, 1.158801, -0.405579,
        -1.278287, 0.391306, 0.857077,

        /* Outliers */
         45, 43, 0,
         23, -3, 2,

        10.6772, 10.7365, 9.9264,
        8.7785, 11.1680, 9.5915,
        8.6872, 9.6464, 10.3801,
        10.0142, 8.8311, 9.2021,
        8.4179, 9.8572, 11.6356,
        9.8487, 10.4979, 10.8620,
        10.0957, 9.7878, 12.2653,
        11.4528, 11.5186, 10.3050,
        10.9284, 9.9654, 10.4562,
        8.5272, 10.7451, 9.8355,
        10.1508, 10.2318, 10.2417,
        10.7342, 10.0689, 9.9918,
        10.4784, 9.2032, 10.6060,
        10.1309, 9.4392, 10.9674,
        10.6971, 10.3347, 11.0447,
        7.9489, 9.4566, 9.5258,
        10.4827, 10.3030, 10.5582,
        10.4496, 10.3880, 11.1661,
        11.0291, 10.0233, 9.9280,
        9.0638, 9.3650, 9.3670,

        /* Outliers */
         32, 54, 23,
         23, 51, 77,


    -34.233,  -30.841,  -31.720,
     -32.629,  -31.786,  -31.290,
     -31.466,  -31.984,  -33.254,
     -31.878,  -33.052,  -31.761,
     -33.528,  -30.921,  -32.836,
     -31.793,  -32.082,  -30.453,
     -31.812,  -32.417,  -31.874,
     -32.127,  -32.599,  -32.806,
     -32.979,  -32.096,  -31.754,
     -31.759,  -31.925,  -31.313,
     -30.531,  -31.838,  -31.179,
     -32.168,  -31.928,  -30.649,
     -31.049,  -32.092,  -31.408,
     -33.006,  -31.753,  -31.961,
     -32.092,  -32.391,  -31.501,
     -31.184,  -31.634,  -32.802,
     -30.658,  -31.616,  -31.493,
     -31.958,  -31.694,  -31.425,
     -33.114,  -32.029,  -31.459,
     -31.081,  -34.486,  -32.020,

        /* Outlier */
        -22, -34, 53 };

    size_t idx[row] = { 0 };
    float epsilon = 10;
    size_t min_pts = 3;

    /* Compute dbscan */
    dbscan(X, idx, epsilon, min_pts, row, column);

    /* Print idx */
    size_t i;
    for (i = 0; i < row; i++) {
        printf("%i\n", idx[i]);
    }

    return 0;
}

输出变为：

1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
0
0
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
0
0
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
3
0
c
Martin Hope
euraad
Asked: 2023-12-28 17:31:41 +0800 CST
  • 5

我试图使conv(x, k, 'same')结果与 相同FFT。但我注意到这仅适用于 shapefull而不适用于same.

% Create vector and kernel
x = randn(1, 100); 
k = randn(1, 10);  

% Do 1D conv
y_conv = conv(x, k, 'full');

% Do 1D FFT
N = length(x) + length(k) - 1;
X = fft(x, N);
K = fft(k, N);
Y_fft = ifft(X .* K);

% Compare
max_diff = max(abs(y_conv - Y_fft));
disp(['The maximal difference is: ' num2str(max_diff)]);

我之所以想要same，是因为我的输出数据需要与输入数据相同。

问题：

为了使 FFT 结果与 相同，需要做什么/更改conv(x, k, 'same')？

matlab
Martin Hope
euraad
Asked: 2023-09-17 06:38:44 +0800 CST
  • 7

假设你有这个矩阵

>> P = zeros(8, 12)
P =

   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0

并且您想添加1每个坐标。如果是同一个坐标，那就是+1。

>> angles = [1 2 1 3 3 2 5 6 2 1 3 8]
angles =

   1   2   1   3   3   2   5   6   2   1   3   8

>> r = [1 2 4 2 5 6 1 4 6 6 3 6]
r =

   1   2   4   2   5   6   1   4   6   6   3   6

如果我使用sub2ind

>> indices = sub2ind(size(P), angles, r);
>> P(indices) = P(indices) + 1;
>> P
P =

   1   0   0   1   0   1   0   0   0   0   0   0
   0   1   0   0   0   1   0   0   0   0   0   0
   0   1   1   0   1   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   1   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   1   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   1   0   0   0   0   0   0

但如果我使用普通索引

>> for k = 1:length(r)
      P(angles(k), r(k)) = P(angles(k), r(k)) + 1;
   end
>> P
P =

   1   0   0   1   0   1   0   0   0   0   0   0
   0   1   0   0   0   2   0   0   0   0   0   0
   0   1   1   0   1   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   1   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   1   0   0   0   0   0   0   0   0
   0   0   0   0   0   0   0   0   0   0   0   0
   0   0   0   0   0   1   0   0   0   0   0   0

仔细观察，2矩阵内部有一个P，但上面的矩阵内部没有P。

问题：

为什么sub2ind索引与常规P(i, j)索引不同？

arrays