Numerical Recipes Forum  

Go Back   Numerical Recipes Forum > Numerical Recipes Third Edition Forum > Methods: All Chapters in NR3

Reply
 
Thread Tools Rating: Thread Rating: 2 votes, 5.00 average. Display Modes
  #1  
Old 04-07-2012, 05:46 PM
mikeme mikeme is offline
Registered User
 
Join Date: Apr 2012
Posts: 2
Chapter 21 - KD Tree

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<3p(x[ii],y[ii],z[ii]);
        
VecPoints[ii] = p;
    }
    
    
// Construct KD Tree Based off of Mean Positions of RSO's
    
KDtree<3tree(VecPoints);
    
tree.nnearest(1,nn,dn,19);
    
cout << &VecPoints[1];
    return 
0;

Reply With Quote
  #2  
Old 04-18-2012, 12:09 PM
mikeme mikeme is offline
Registered User
 
Join Date: Apr 2012
Posts: 2
New Question

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?
Reply With Quote
  #3  
Old 04-21-2012, 01:16 AM
davekw7x davekw7x is offline
Registered User
 
Join Date: Jan 2008
Posts: 453
Quote:
Originally Posted by mikeme View Post
.... 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.


Code:
#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:

Code:
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

Last edited by davekw7x; 04-21-2012 at 08:53 AM.
Reply With Quote
  #4  
Old 02-02-2013, 06:32 PM
amirvahid amirvahid is offline
Registered User
 
Join Date: Mar 2008
Location: Notre Dame IN USA
Posts: 67
Data file

Could you guys upload your data files so everyone can benefit from your post? Thanks!
Reply With Quote
  #5  
Old 02-03-2013, 07:16 AM
davekw7x davekw7x is offline
Registered User
 
Join Date: Jan 2008
Posts: 453
Quote:
Originally Posted by amirvahid View Post
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:
Code:
//
// 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
Reply With Quote
Reply

Tags
chapter 21, computational geometry, kd tree, nearest neighbors

Thread Tools
Display Modes Rate This Thread
Rate This Thread:

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off

Forum Jump


All times are GMT -5. The time now is 11:14 PM.


Powered by vBulletin® Version 3.8.6
Copyright ©2000 - 2014, Jelsoft Enterprises Ltd.