C# kdtree类代码示例

        /*************************************************************************
        KD-tree creation

        This subroutine creates KD-tree from set of X-values and optional Y-values

        INPUT PARAMETERS
            XY      -   dataset, array[0..N-1,0..NX+NY-1].
                        one row corresponds to one point.
                        first NX columns contain X-values, next NY (NY may be zero)
                        columns may contain associated Y-values
            N       -   number of points, N>=1
            NX      -   space dimension, NX>=1.
            NY      -   number of optional Y-values, NY>=0.
            NormType-   norm type:
                        * 0 denotes infinity-norm
                        * 1 denotes 1-norm
                        * 2 denotes 2-norm (Euclidean norm)
                        
        OUTPUT PARAMETERS
            KDT     -   KD-tree
            
            
        NOTES

        1. KD-tree  creation  have O(N*logN) complexity and O(N*(2*NX+NY))  memory
           requirements.
        2. Although KD-trees may be used with any combination of N  and  NX,  they
           are more efficient than brute-force search only when N >> 4^NX. So they
           are most useful in low-dimensional tasks (NX=2, NX=3). NX=1  is another
           inefficient case, because  simple  binary  search  (without  additional
           structures) is much more efficient in such tasks than KD-trees.

          -- ALGLIB --
             Copyright 28.02.2010 by Bochkanov Sergey
        *************************************************************************/
        public static void kdtreebuild(ref double[,] xy,
            int n,
            int nx,
            int ny,
            int normtype,
            ref kdtree kdt)
        {
            int[] tags = new int[0];
            int i = 0;

            System.Diagnostics.Debug.Assert(n>=1, "KDTreeBuild: N<1!");
            System.Diagnostics.Debug.Assert(nx>=1, "KDTreeBuild: NX<1!");
            System.Diagnostics.Debug.Assert(ny>=0, "KDTreeBuild: NY<0!");
            System.Diagnostics.Debug.Assert(normtype>=0 & normtype<=2, "KDTreeBuild: incorrect NormType!");
            tags = new int[n];
            for(i=0; i<=n-1; i++)
            {
                tags[i] = 0;
            }
            kdtreebuildtagged(ref xy, ref tags, n, nx, ny, normtype, ref kdt);
        }

        /*************************************************************************
        Copies X[] to KDT.X[]
        Loads distance from X[] to bounding box.
        Initializes CurBox[].

          -- ALGLIB --
             Copyright 28.02.2010 by Bochkanov Sergey
        *************************************************************************/
        private static void kdtreeinitbox(kdtree kdt,
            double[] x)
        {
            int i = 0;
            double vx = 0;
            double vmin = 0;
            double vmax = 0;

            alglib.ap.assert(kdt.n>0, "KDTreeInitBox: internal error");
            
            //
            // calculate distance from point to current bounding box
            //
            kdt.curdist = 0;
            if( kdt.normtype==0 )
            {
                for(i=0; i<=kdt.nx-1; i++)
                {
                    vx = x[i];
                    vmin = kdt.boxmin[i];
                    vmax = kdt.boxmax[i];
                    kdt.x[i] = vx;
                    kdt.curboxmin[i] = vmin;
                    kdt.curboxmax[i] = vmax;
                    if( (double)(vx)<(double)(vmin) )
                    {
                        kdt.curdist = Math.Max(kdt.curdist, vmin-vx);
                    }
                    else
                    {
                        if( (double)(vx)>(double)(vmax) )
                        {
                            kdt.curdist = Math.Max(kdt.curdist, vx-vmax);
                        }
                    }
                }
            }
            if( kdt.normtype==1 )
            {
                for(i=0; i<=kdt.nx-1; i++)
                {
                    vx = x[i];
                    vmin = kdt.boxmin[i];
                    vmax = kdt.boxmax[i];
                    kdt.x[i] = vx;
                    kdt.curboxmin[i] = vmin;
                    kdt.curboxmax[i] = vmax;
                    if( (double)(vx)<(double)(vmin) )
                    {
                        kdt.curdist = kdt.curdist+vmin-vx;
                    }
                    else
                    {
                        if( (double)(vx)>(double)(vmax) )
                        {
                            kdt.curdist = kdt.curdist+vx-vmax;
                        }
                    }
                }
            }
            if( kdt.normtype==2 )
            {
                for(i=0; i<=kdt.nx-1; i++)
                {
                    vx = x[i];
                    vmin = kdt.boxmin[i];
                    vmax = kdt.boxmax[i];
                    kdt.x[i] = vx;
                    kdt.curboxmin[i] = vmin;
                    kdt.curboxmax[i] = vmax;
                    if( (double)(vx)<(double)(vmin) )
                    {
                        kdt.curdist = kdt.curdist+math.sqr(vmin-vx);
                    }
                    else
                    {
                        if( (double)(vx)>(double)(vmax) )
                        {
                            kdt.curdist = kdt.curdist+math.sqr(vx-vmax);
                        }
                    }
                }
            }
        }

        /*************************************************************************
        This function allocates all dataset-dependent array fields of KDTree, i.e.
        such array fields that their dimensions depend on dataset size.

        This function do not sets KDT.N, KDT.NX or KDT.NY -
        it just allocates arrays.

          -- ALGLIB --
             Copyright 14.03.2011 by Bochkanov Sergey
        *************************************************************************/
        private static void kdtreeallocdatasetdependent(kdtree kdt,
            int n,
            int nx,
            int ny)
        {
            alglib.ap.assert(n>0, "KDTreeAllocDatasetDependent: internal error");
            kdt.xy = new double[n, 2*nx+ny];
            kdt.tags = new int[n];
            kdt.idx = new int[n];
            kdt.r = new double[n];
            kdt.x = new double[nx];
            kdt.buf = new double[Math.Max(n, nx)];
            kdt.nodes = new int[splitnodesize*2*n];
            kdt.splits = new double[2*n];
        }

        /*************************************************************************
        Serializer: unserialization

          -- ALGLIB --
             Copyright 14.03.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void kdtreeunserialize(alglib.serializer s,
            kdtree tree)
        {
            int i0 = 0;
            int i1 = 0;

            
            //
            // check correctness of header
            //
            i0 = s.unserialize_int();
            alglib.ap.assert(i0==scodes.getkdtreeserializationcode(), "KDTreeUnserialize: stream header corrupted");
            i1 = s.unserialize_int();
            alglib.ap.assert(i1==kdtreefirstversion, "KDTreeUnserialize: stream header corrupted");
            
            //
            // Unserialize data
            //
            tree.n = s.unserialize_int();
            tree.nx = s.unserialize_int();
            tree.ny = s.unserialize_int();
            tree.normtype = s.unserialize_int();
            apserv.unserializerealmatrix(s, ref tree.xy);
            apserv.unserializeintegerarray(s, ref tree.tags);
            apserv.unserializerealarray(s, ref tree.boxmin);
            apserv.unserializerealarray(s, ref tree.boxmax);
            apserv.unserializeintegerarray(s, ref tree.nodes);
            apserv.unserializerealarray(s, ref tree.splits);
            kdtreealloctemporaries(tree, tree.n, tree.nx, tree.ny);
        }

        /*************************************************************************
        Recursive kd-tree generation subroutine.

        PARAMETERS
            KDT         tree
            NodesOffs   unused part of Nodes[] which must be filled by tree
            SplitsOffs  unused part of Splits[]
            I1, I2      points from [I1,I2) are processed
            
        NodesOffs[] and SplitsOffs[] must be large enough.

          -- ALGLIB --
             Copyright 28.02.2010 by Bochkanov Sergey
        *************************************************************************/
        private static void kdtreegeneratetreerec(kdtree kdt,
            ref int nodesoffs,
            ref int splitsoffs,
            int i1,
            int i2,
            int maxleafsize)
        {
            int n = 0;
            int nx = 0;
            int ny = 0;
            int i = 0;
            int j = 0;
            int oldoffs = 0;
            int i3 = 0;
            int cntless = 0;
            int cntgreater = 0;
            double minv = 0;
            double maxv = 0;
            int minidx = 0;
            int maxidx = 0;
            int d = 0;
            double ds = 0;
            double s = 0;
            double v = 0;
            int i_ = 0;
            int i1_ = 0;

            alglib.ap.assert(kdt.n>0, "KDTreeGenerateTreeRec: internal error");
            alglib.ap.assert(i2>i1, "KDTreeGenerateTreeRec: internal error");
            
            //
            // Generate leaf if needed
            //
            if( i2-i1<=maxleafsize )
            {
                kdt.nodes[nodesoffs+0] = i2-i1;
                kdt.nodes[nodesoffs+1] = i1;
                nodesoffs = nodesoffs+2;
                return;
            }
            
            //
            // Load values for easier access
            //
            nx = kdt.nx;
            ny = kdt.ny;
            
            //
            // select dimension to split:
            // * D is a dimension number
            //
            d = 0;
            ds = kdt.curboxmax[0]-kdt.curboxmin[0];
            for(i=1; i<=nx-1; i++)
            {
                v = kdt.curboxmax[i]-kdt.curboxmin[i];
                if( (double)(v)>(double)(ds) )
                {
                    ds = v;
                    d = i;
                }
            }
            
            //
            // Select split position S using sliding midpoint rule,
            // rearrange points into [I1,I3) and [I3,I2)
            //
            s = kdt.curboxmin[d]+0.5*ds;
            i1_ = (i1) - (0);
            for(i_=0; i_<=i2-i1-1;i_++)
            {
                kdt.buf[i_] = kdt.xy[i_+i1_,d];
            }
            n = i2-i1;
            cntless = 0;
            cntgreater = 0;
            minv = kdt.buf[0];
            maxv = kdt.buf[0];
            minidx = i1;
            maxidx = i1;
            for(i=0; i<=n-1; i++)
            {
                v = kdt.buf[i];
                if( (double)(v)<(double)(minv) )
                {
                    minv = v;
//.........这里部分代码省略.........

        /*************************************************************************
        Tags  from  last  query;  'interactive' variant for languages like  Python
        which  support  constructs  like "Tags = KDTreeQueryResultsTagsI(KDT)" and
        interactive mode of interpreter.

        This function allocates new array on each call,  so  it  is  significantly
        slower than its 'non-interactive' counterpart, but it is  more  convenient
        when you call it from command line.

          -- ALGLIB --
             Copyright 28.02.2010 by Bochkanov Sergey
        *************************************************************************/
        public static void kdtreequeryresultstagsi(kdtree kdt,
            ref int[] tags)
        {
            tags = new int[0];

            kdtreequeryresultstags(kdt, ref tags);
        }

        /*************************************************************************
        Serializer: allocation

          -- ALGLIB --
             Copyright 14.03.2011 by Bochkanov Sergey
        *************************************************************************/
        public static void kdtreealloc(alglib.serializer s,
            kdtree tree)
        {
            
            //
            // Header
            //
            s.alloc_entry();
            s.alloc_entry();
            
            //
            // Data
            //
            s.alloc_entry();
            s.alloc_entry();
            s.alloc_entry();
            s.alloc_entry();
            apserv.allocrealmatrix(s, tree.xy, -1, -1);
            apserv.allocintegerarray(s, tree.tags, -1);
            apserv.allocrealarray(s, tree.boxmin, -1);
            apserv.allocrealarray(s, tree.boxmax, -1);
            apserv.allocintegerarray(s, tree.nodes, -1);
            apserv.allocrealarray(s, tree.splits, -1);
        }

            /*************************************************************************
            This function allocates temporaries.

            This function do not sets KDT.N, KDT.NX or KDT.NY -
            it just allocates arrays.

              -- ALGLIB --
                 Copyright 14.03.2011 by Bochkanov Sergey
            *************************************************************************/
            private static void kdtreealloctemporaries(kdtree kdt,
                int n,
                int nx,
                int ny) {
                kdt.x = new double[nx];
                kdt.idx = new int[n];
                kdt.r = new double[n];
                kdt.buf = new double[Math.Max(n, nx)];
                kdt.curboxmin = new double[nx];
                kdt.curboxmax = new double[nx];
            }

 public override alglib.apobject make_copy()
 {
     kdtree _result = new kdtree();
     _result.n = n;
     _result.nx = nx;
     _result.ny = ny;
     _result.normtype = normtype;
     _result.xy = (double[,])xy.Clone();
     _result.tags = (int[])tags.Clone();
     _result.boxmin = (double[])boxmin.Clone();
     _result.boxmax = (double[])boxmax.Clone();
     _result.nodes = (int[])nodes.Clone();
     _result.splits = (double[])splits.Clone();
     _result.x = (double[])x.Clone();
     _result.kneeded = kneeded;
     _result.rneeded = rneeded;
     _result.selfmatch = selfmatch;
     _result.approxf = approxf;
     _result.kcur = kcur;
     _result.idx = (int[])idx.Clone();
     _result.r = (double[])r.Clone();
     _result.buf = (double[])buf.Clone();
     _result.curboxmin = (double[])curboxmin.Clone();
     _result.curboxmax = (double[])curboxmax.Clone();
     _result.curdist = curdist;
     _result.debugcounter = debugcounter;
     return _result;
 }

            /*************************************************************************
            This function allocates all dataset-independent array  fields  of  KDTree,
            i.e.  such  array  fields  that  their dimensions do not depend on dataset
            size.

            This function do not sets KDT.NX or KDT.NY - it just allocates arrays

              -- ALGLIB --
                 Copyright 14.03.2011 by Bochkanov Sergey
            *************************************************************************/
            private static void kdtreeallocdatasetindependent(kdtree kdt,
                int nx,
                int ny) {
                kdt.x = new double[nx];
                kdt.boxmin = new double[nx];
                kdt.boxmax = new double[nx];
                kdt.curboxmin = new double[nx];
                kdt.curboxmax = new double[nx];
            }

            /*************************************************************************
            This function allocates all dataset-dependent array fields of KDTree, i.e.
            such array fields that their dimensions depend on dataset size.

            This function do not sets KDT.N, KDT.NX or KDT.NY -
            it just allocates arrays.

              -- ALGLIB --
                 Copyright 14.03.2011 by Bochkanov Sergey
            *************************************************************************/
            private static void kdtreeallocdatasetdependent(kdtree kdt,
                int n,
                int nx,
                int ny) {
                kdt.xy = new double[n, 2 * nx + ny];
                kdt.tags = new int[n];
                kdt.idx = new int[n];
                kdt.r = new double[n];
                kdt.x = new double[nx];
                kdt.buf = new double[Math.Max(n, nx)];
                kdt.nodes = new int[splitnodesize * 2 * n];
                kdt.splits = new double[2 * n];
            }

            /*************************************************************************
            X-values from last query

            INPUT PARAMETERS
                KDT     -   KD-tree
                X       -   possibly pre-allocated buffer. If X is too small to store
                            result, it is resized. If size(X) is enough to store
                            result, it is left unchanged.

            OUTPUT PARAMETERS
                X       -   rows are filled with X-values

            NOTES
            1. points are ordered by distance from the query point (first = closest)
            2. if  XY is larger than required to store result, only leading part  will
               be overwritten; trailing part will be left unchanged. So  if  on  input
               XY = [[A,B],[C,D]], and result is [1,2],  then  on  exit  we  will  get
               XY = [[1,2],[C,D]]. This is done purposely to increase performance;  if
               you want function  to  resize  array  according  to  result  size,  use
               function with same name and suffix 'I'.

            SEE ALSO
            * KDTreeQueryResultsXY()            X- and Y-values
            * KDTreeQueryResultsTags()          tag values
            * KDTreeQueryResultsDistances()     distances

              -- ALGLIB --
                 Copyright 28.02.2010 by Bochkanov Sergey
            *************************************************************************/
            public static void kdtreequeryresultsx(kdtree kdt,
                ref double[,] x) {
                int i = 0;
                int k = 0;
                int i_ = 0;
                int i1_ = 0;

                if(kdt.kcur == 0) {
                    return;
                }
                if(ap.rows(x) < kdt.kcur | ap.cols(x) < kdt.nx) {
                    x = new double[kdt.kcur, kdt.nx];
                }
                k = kdt.kcur;
                for(i = 0; i <= k - 1; i++) {
                    i1_ = (kdt.nx) - (0);
                    for(i_ = 0; i_ <= kdt.nx - 1; i_++) {
                        x[i, i_] = kdt.xy[kdt.idx[i], i_ + i1_];
                    }
                }
            }

            /*************************************************************************
            R-NN query: all points within R-sphere centered at X

            INPUT PARAMETERS
                KDT         -   KD-tree
                X           -   point, array[0..NX-1].
                R           -   radius of sphere (in corresponding norm), R>0
                SelfMatch   -   whether self-matches are allowed:
                                * if True, nearest neighbor may be the point itself
                                  (if it exists in original dataset)
                                * if False, then only points with non-zero distance
                                  are returned
                                * if not given, considered True

            RESULT
                number of neighbors found, >=0

            This  subroutine  performs  query  and  stores  its result in the internal
            structures of the KD-tree. You can use  following  subroutines  to  obtain
            actual results:
            * KDTreeQueryResultsX() to get X-values
            * KDTreeQueryResultsXY() to get X- and Y-values
            * KDTreeQueryResultsTags() to get tag values
            * KDTreeQueryResultsDistances() to get distances

              -- ALGLIB --
                 Copyright 28.02.2010 by Bochkanov Sergey
            *************************************************************************/
            public static int kdtreequeryrnn(kdtree kdt,
                double[] x,
                double r,
                bool selfmatch) {
                int result = 0;
                int i = 0;
                int j = 0;

                ap.assert((double)(r) > (double)(0), "KDTreeQueryRNN: incorrect R!");
                ap.assert(ap.len(x) >= kdt.nx, "KDTreeQueryRNN: Length(X)<NX!");
                ap.assert(apserv.isfinitevector(x, kdt.nx), "KDTreeQueryRNN: X contains infinite or NaN values!");

                //
                // Prepare parameters
                //
                kdt.kneeded = 0;
                if(kdt.normtype != 2) {
                    kdt.rneeded = r;
                } else {
                    kdt.rneeded = math.sqr(r);
                }
                kdt.selfmatch = selfmatch;
                kdt.approxf = 1;
                kdt.kcur = 0;

                //
                // calculate distance from point to current bounding box
                //
                kdtreeinitbox(kdt, x);

                //
                // call recursive search
                // results are returned as heap
                //
                kdtreequerynnrec(kdt, 0);

                //
                // pop from heap to generate ordered representation
                //
                // last element is not pop'ed because it is already in
                // its place
                //
                result = kdt.kcur;
                j = kdt.kcur;
                for(i = kdt.kcur; i >= 2; i--) {
                    tsort.tagheappopi(ref kdt.r, ref kdt.idx, ref j);
                }
                return result;
            }

        /*************************************************************************
        KD-tree creation

        This  subroutine  creates  KD-tree  from set of X-values, integer tags and
        optional Y-values

        INPUT PARAMETERS
            XY      -   dataset, array[0..N-1,0..NX+NY-1].
                        one row corresponds to one point.
                        first NX columns contain X-values, next NY (NY may be zero)
                        columns may contain associated Y-values
            Tags    -   tags, array[0..N-1], contains integer tags associated
                        with points.
            N       -   number of points, N>=1
            NX      -   space dimension, NX>=1.
            NY      -   number of optional Y-values, NY>=0.
            NormType-   norm type:
                        * 0 denotes infinity-norm
                        * 1 denotes 1-norm
                        * 2 denotes 2-norm (Euclidean norm)

        OUTPUT PARAMETERS
            KDT     -   KD-tree

        NOTES

        1. KD-tree  creation  have O(N*logN) complexity and O(N*(2*NX+NY))  memory
           requirements.
        2. Although KD-trees may be used with any combination of N  and  NX,  they
           are more efficient than brute-force search only when N >> 4^NX. So they
           are most useful in low-dimensional tasks (NX=2, NX=3). NX=1  is another
           inefficient case, because  simple  binary  search  (without  additional
           structures) is much more efficient in such tasks than KD-trees.

          -- ALGLIB --
             Copyright 28.02.2010 by Bochkanov Sergey
        *************************************************************************/
        public static void kdtreebuildtagged(double[,] xy,
            int[] tags,
            int n,
            int nx,
            int ny,
            int normtype,
            kdtree kdt)
        {
            int i = 0;
            int j = 0;
            int maxnodes = 0;
            int nodesoffs = 0;
            int splitsoffs = 0;
            int i_ = 0;
            int i1_ = 0;

            ap.assert(n>=1, "KDTreeBuildTagged: N<1!");
            ap.assert(nx>=1, "KDTreeBuildTagged: NX<1!");
            ap.assert(ny>=0, "KDTreeBuildTagged: NY<0!");
            ap.assert(normtype>=0 & normtype<=2, "KDTreeBuildTagged: incorrect NormType!");
            ap.assert(ap.rows(xy)>=n, "KDTreeBuildTagged: rows(X)<N!");
            ap.assert(ap.cols(xy)>=nx+ny, "KDTreeBuildTagged: cols(X)<NX+NY!");
            ap.assert(apserv.apservisfinitematrix(xy, n, nx+ny), "KDTreeBuildTagged: X contains infinite or NaN values!");
            
            //
            // initialize
            //
            kdt.n = n;
            kdt.nx = nx;
            kdt.ny = ny;
            kdt.normtype = normtype;
            kdt.distmatrixtype = 0;
            kdt.xy = new double[n, 2*nx+ny];
            kdt.tags = new int[n];
            kdt.idx = new int[n];
            kdt.r = new double[n];
            kdt.x = new double[nx];
            kdt.buf = new double[Math.Max(n, nx)];
            
            //
            // Initial fill
            //
            for(i=0; i<=n-1; i++)
            {
                for(i_=0; i_<=nx-1;i_++)
                {
                    kdt.xy[i,i_] = xy[i,i_];
                }
                i1_ = (0) - (nx);
                for(i_=nx; i_<=2*nx+ny-1;i_++)
                {
                    kdt.xy[i,i_] = xy[i,i_+i1_];
                }
                kdt.tags[i] = tags[i];
            }
            
            //
            // Determine bounding box
            //
            kdt.boxmin = new double[nx];
            kdt.boxmax = new double[nx];
            kdt.curboxmin = new double[nx];
            kdt.curboxmax = new double[nx];
//.........这里部分代码省略.........

        /*************************************************************************
        X-values from last query; 'interactive' variant for languages like  Python
        which   support    constructs   like  "X = KDTreeQueryResultsXI(KDT)"  and
        interactive mode of interpreter.

        This function allocates new array on each call,  so  it  is  significantly
        slower than its 'non-interactive' counterpart, but it is  more  convenient
        when you call it from command line.

          -- ALGLIB --
             Copyright 28.02.2010 by Bochkanov Sergey
        *************************************************************************/
        public static void kdtreequeryresultsxi(kdtree kdt,
            ref double[,] x)
        {
            x = new double[0,0];

            kdtreequeryresultsx(kdt, ref x);
        }