Chapter 21 - KD Tree

[mikeme]

Hi all,

I'm very new to C++, but need to implement a kd tree for some research. I believe I've implemented the code correctly, however, I cannot obtain the output of the result of the nearest neighbor search. Additionally, in regards to the external file I've loaded in, does the selecti() routine alter the order of the points? If so, after the nearest neighbor indices are obtained, how do I match the actual point locations to each index? Any help would be much appreciated. I'm trying very hard to understand this, but I'm limited on time.

PHP Code:

#include <iostream>
#include <fstream>
#include <sstream>
#include <string>
#include <vector>
#include "nr3.h"
#include "pointbox.h"
#include "kdtree.h"
using namespace std;

Int main(){
    // Initialize Variables
    vector < Point<3> > VecPoints(2000);
    Doub x[2000];
    Doub y[2000];
    Doub z[2000];
    Int i=0;
    Int j=0;
    Int k=0;
    Int *nn =new Int;
    Doub *dn=new Doub;
    
    // Load x-position data
    ifstream infile1("meanposx1.txt", ios::in);
    while(! infile1.eof())
    {
        infile1 >> x[i];
        i++;
    }
    infile1.close();
    
    // Load y-position data
    ifstream infile2("meanposy1.txt", ios::in);
    while(! infile2.eof())
    {
        infile2 >> y[j];
        j++;
    }
    infile2.close();
    
    // Load z-position data
    ifstream infile3("meanposz1.txt", ios::in);
    while(! infile3.eof())
    {
        infile3 >> z[k];
        k++;
    }
    infile3.close();
    
    // Create Vector of Points
    for (Int ii=0;ii<2000;ii++)
    {
        Point<3> p(x[ii],y[ii],z[ii]);
        VecPoints[ii] = p;
    }
    
    // Construct KD Tree Based off of Mean Positions of RSO's
    KDtree<3> tree(VecPoints);
    tree.nnearest(1,nn,dn,19);
    cout << &VecPoints[1];
    return 0;
} 


[mikeme]

Ok, so I've made my way through creating the tree and finding the nearest neighbors. My question is that when I ask for, say, 5 nearest neighbors, the search only returns the fifth nearest neighbor as opposed to returning a list of the five closest. Can anyone help me out with this portion?

[davekw7x]

Originally Posted by mikeme

.... Can anyone help me out ...

You must create storage for arrays that can hold the index values and distances for the number of neighbors that you tell it to compute.

#include "../code/nr3.h"
#include "../code/pointbox.h"
#include "../code/kdtree.h"
//
// A convenient function to print coordinates of any kind of NR Point.
//
template <int DIM>
void printCoords(const Point<DIM> & pt);

int main(int argc, char **argv)
{
    //
    // Your stuff to read the coordinates from the file and
    // create vecPoints, a vector of Point<3> objects.
    // Make sure it reads all of the values that you expect.
    //
    // Then...

    int numPoints = vecPoints.size();
    cout << fixed;
    cout << "Number of points read = " << numPoints << endl;

    // Create a KDtree object 
    KDtree <3> tree(vecPoints);

    //
    // Get the index value of the point that you are investigating.
    // If you enter a non-numeric value here, the program
    // terminates
    //
    Int indx;
    cout << "Enter the index of a point in the tree: ";
    while (cin >> indx) {
        if (indx < 0 || indx >= numPoints){
            cout << "Must be an integer 0..." << numPoints-1 << endl;
        }
        else {
            // Find nearest five points
            const int numNearest = 5;
            Int nn[numNearest];
            Doub dn[numNearest];

            tree.nnearest(indx, nn, dn, numNearest);

            cout << endl << "Reference point "
                 << setw(4) << indx << ": ";
            printCoords(vecPoints[indx]);
            cout << endl;

            cout << "Here are the " << numNearest
                 << " points that are nearest to point "
                 << indx << endl << endl;

            cout << "Index  Distance     Coordinates" << endl;
            cout << "-----------------------------------------------" << endl;
            
            for (int i = 0; i < numNearest; i++) {
                cout << setw(4) << nn[i] << "   " << dn[i]
                     << "    ";
                printCoords(vecPoints[nn[i]]);
                cout << endl;
            }
        }
        // Go back and try again with another point
        cout << endl;
        cout << "Enter the index of another input point: ";
    }
    cout << "Goodbye for now." << endl;

    return 0;
}

template <int DIM>
void printCoords(const Point<DIM> & pt)
{
    cout << "(";
    for (int i = 0; i < DIM; i++) {
        cout << pt.x[i];
        if (i < DIM-1) {
            cout << ",";
        }
    }
    cout << ")";
}

My file had points that were inside a cube with coordinates between 0 and 1. Here's the output for a run.

Program output is blue; user input is black:

Number of points read = 2000
Enter the index of an input point: 0

Reference point    0: (0.035093,0.767641,0.628737)
Here are the 5 points that are nearest to point 0

Index  Distance     Coordinates
-----------------------------------------------
 803   0.085480    (0.028558,0.807263,0.553278)
 836   0.075850    (0.067090,0.834319,0.645573)
1939   0.074550    (0.108575,0.778370,0.622174)
 202   0.066690    (0.065788,0.723138,0.589687)
 293   0.074821    (0.019028,0.840708,0.627655)

Enter the index of another input point: 1999

Reference point 1999: (0.885194,0.391274,0.464631)
Here are the 5 points that are nearest to point 1999

Index  Distance     Coordinates
-----------------------------------------------
1083   0.075079    (0.817496,0.361783,0.451058)
1234   0.074557    (0.824474,0.374998,0.504718)
1736   0.044114    (0.898018,0.349300,0.460174)
1670   0.026369    (0.863950,0.404000,0.473690)
1325   0.045187    (0.885863,0.346117,0.463137)

Enter the index of another input point: quit
Goodbye for now.

Regards,

Dave

[amirvahid]

Could you guys upload your data files so everyone can benefit from your post? Thanks!

[davekw7x]

Originally Posted by amirvahid

Could you guys upload your data files so everyone can benefit from your post? Thanks!

Here's how I generated the file that gave the output from my previous post:

//
// genboxpoints.cpp
//
// Generates 2000 3-D points.
//
// Points are inside a unit cube.
//
//  davekw7x
//
#include "../code/nr3.h"
#include "../code/ran.h"

int main()
{
    const char * out_name = "2000points.dat";
    ofstream out_file(out_name);
    if (!out_file) {
        cerr << "Can't open file " << out_name << " for writing." << endl;
        exit(EXIT_FAILURE);
    }
    cout << "Writing to file " << out_name << endl;

    Ran ran(123456);
    out_file << scientific;

    for (int i = 0; i < 2000; i++) {
        out_file << ran.doub() << "  " << ran.doub() << "  " << ran.doub() << endl;
    }

    out_file.close();

    return 0;
}

Regards,

Dave