Найти неэффективных ближайших соседей в трехмерном пространстве

Question

Мой входной файл представляет собой CSV-файл со структурой данных следующим образом:

5,3.5644,5.4556,3.5665
...
int_id,x_float,y_float,z_float

У меня есть структура, которая содержит идентификатор точки и ее три координаты.Мне нужно найти 4 ближайших структуры на основе евклидова расстояния.Я сделал это наивным подходом, но есть ли эффективный способ реализовать это?Я читал об алгоритме knn, но для этого нужны внешние библиотеки.

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <assert.h>

//struct for input csv data
struct oxygen_coordinates
{
    unsigned int index; //index of an atom
    //x,y and z coordinates of atom 
    float x;
    float y;
    float z;
};

//Given the data in a line in a an inputfile, process it to put it in   oxygen_coordinates struct
struct oxygen_coordinates * line_splitter(struct oxygen_coordinates *data, const char *input)
{
    return (sscanf(input, "%u,%f,%f,%f", &data->index, &data->x, &data->y, &data->z) != 7)
            ? NULL : data;
}
//Distance function for two pints in a struct
float getDistance(struct oxygen_coordinates a, struct oxygen_coordinates b)
{
      float distance;
      distance = sqrt((a.x - b.x) * (a.x - b.x) + (a.y-b.y) *(a.y-b.y) + (a.z - b.z) * (a.z - b.z));
      return distance;
} 

//struct for neighbour distance and their indices
    struct nbrs_ind {
    float value;
    int index;
    };

// comparision function for sorting the neighbours -> qsort library
int cmp(const void *pa, const void *pb)
{
    struct nbrs_ind *pa1 = (struct nbrs_ind *)pa;
    struct nbrs_ind *pa2 = (struct nbrs_ind *)pb;
    if ((*pa1).value < (*pa2).value)
            return -1;
    else if ((*pa1).value > (*pa2).value)
            return 1;
    else
            return 0;
}
//main program
int  main(int argc, char *argv[])
{
    FILE *stream; // file pointer
    char *line = NULL; //line pointer
    size_t len = 0;
    ssize_t nread;
    struct oxygen_coordinates * atom_data = NULL; //pointer to oxygen_coordinate struct
    int numatoms = 0; // counter variable for number of atoms

    int i,j,k,p ; // loop initilizers
    //Check for correct number of arguments
    if (argc !=2 ) {
            fprintf(stderr, "Usage: %s <inputfile> <outputfile>\n", argv[0]);
            exit(EXIT_FAILURE);
    }
    //Open the input csv file given in the first argument
    stream = fopen(argv[1], "r");
    if (stream == NULL) {
            perror("fopen");
            exit(EXIT_FAILURE);
    }
    while ((nread = getline(&line, &len, stream)) != -1) {
            if ((atom_data = realloc(atom_data, (size_t) (numatoms + 1) * sizeof(struct oxygen_coordinates))) == NULL) {
                    fprintf(stderr, "error not enough memory");
                    exit(EXIT_FAILURE);
            }
            line_splitter(&atom_data[numatoms], line);
            numatoms = numatoms + 1;
    }
    free(line);
    fclose(stream);

    // All the data is read in memory in atom_data
    printf("-------------------------------------------\n");
    printf("There are %d atoms in the input file. \n", numatoms);
    printf("-------------------------------------------\n");


   // declare a global array that will hold the 4 nearest atom_data...
   float dist_mat[numatoms][numatoms] ;// create n by n matrix for  distances
    // Create a 2-D distnace matrix        
    for(j=0; j < numatoms; j++){
            for(k=0; k < numatoms; k++) {
                    dist_mat[j][k] = getDistance(atom_data[j], atom_data[k]);
                    printf("%f\t", dist_mat[j][k]);
            }
            printf("\n");
    }
    //now I sort every row from dist_mat and get the closest 4 
    // I need something like as follows
    ////knn(atom_data[query],atom_data,4);//this should return closest 4 points based on  Euclidean distances in atom_data
    free(atom_data);
    exit(EXIT_SUCCESS);
}

Answer 1

Один из способов улучшить производительность - это понять, что вам не нужно фактическое расстояние.Сравнение квадрата расстояния достаточно хорошо, поэтому вы можете пропустить вызов функции sqrt.

Еще одна вещь, которая может, но не обязательно, ускорить процесс, - начать с вычисления только расстояния x.Используйте тот факт, что расстояния всегда положительны, поэтому, если расстояние х длиннее, чем общее расстояние четвертой ближайшей точки, вычислять (a.y-b.y) *(a.y-b.y) + (a.z - b.z) * (a.z - b.z).

не нужно. Если вы выбираете подход, начинающийся только сх-значения, я бы также предложил изменить структуры данных.Вместо структуры для каждой точки используйте четыре массива: int *x, *y, *z, *indexes; Это сделает код более удобным для кеша. (И да, между указателями и массивами есть различия, но здесь это было не так важно)

Приведенные выше методы довольно просты.Если вы хотите продвинуться дальше, вы можете использовать эту идею.

1. Разделите пространство на сетку 3D-блоков.
2. Вычислите, какие точки принадлежат к данному блоку, и сохранитеинформация.

Посмотрите на это изображение:

Например, если у вас есть точка в D4, и выЕсли вы хотите найти ближайших четырех соседей, и вы найдете соседа в D4, то вы знаете , что за пределами квадрата C3: E5 не может быть более близкого соседа.Точно так же точка в D4, у которой есть сосед в D3, не может иметь более близкого соседа за пределами области B3: F6.

Но первое, что нужно при оптимизации, это ВСЕГДА определить узкое место.Вы уверены, что это проблема?Вы говорите, что читаете данные из файла, и чтение одной строки из файла ДОЛЖНО быть медленнее, чем вычисление расстояния.

Answer 2

Помимо упомянутой неэффективности квадратного корня, эти циклы вычисляют расстояние каждой пары дважды .

for(j=0; j < numatoms; j++){
    for(k=0; k < numatoms; k++) {
        dist_mat[j][k] = getDistance(atom_data[j], atom_data[k]);
        printf("%f\t", dist_mat[j][k]);
    }
    printf("\n");
}

Вы можете уменьшить вдвое время выполнения, вычислив расстояние каждой пары только один раз , как

for(j=0; j < numatoms-1; j++){
    for(k=j+1; k < numatoms; k++) {
        dist_mat[j][k] = dist_mat[k][j] = getDistance(atom_data[j], atom_data[k]);
        printf("%f\t", dist_mat[j][k]);
    }
    printf("\n");
}

И, конечно,, когда j == k расстояние должно быть 0.