Author Topic: {Challenge} connect each two points whose distance is shorter than given value  (Read 100648 times)

0 Members and 1 Guest are viewing this topic.

It's Alive!

  • Retired
  • Needs a day job
  • Posts: 8662
  • AKA Daniel
brute force with arx

doit
GetDist: 70
1.036883 Seconds

Code: [Select]
static void ArxDist_doit(void)
  {
    double distance = 0.0;
    acedGetDist(0,ACRX_T("\nGetDist: "),&distance);

    __int64 freq, start, end, diff;
    QueryPerformanceFrequency((LARGE_INTEGER*)&freq);
    QueryPerformanceCounter((LARGE_INTEGER*)&start);

    std::vector<AcGePoint3d> points;
    AcDbBlockTableRecord *pRec;
    AcDbDatabase *pDb = acdbHostApplicationServices()->workingDatabase();
   
    if(acdbOpenAcDbObject((AcDbObject *&) pRec,pDb->currentSpaceId(), AcDb::kForWrite) == Acad::eOk)
    {
      AcDbBlockTableRecordIterator *pIter = NULL;
      if(pRec->newIterator(pIter) == Acad::eOk)
      {
        for(pIter->start(); !pIter->done(); pIter->step())
        {
          AcDbEntity *pEnt = NULL;
          if(pIter->getEntity(pEnt, AcDb::kForRead) == Acad::eOk)
          {
            AcDbPoint *pPoint = AcDbPoint::cast(pEnt);
            if(pPoint)
              points.push_back(pPoint->position());
            pEnt->close();
          }
        }
        delete pIter;
      }
      size_t len = points.size();
      for (size_t i = 0; i < len; i++)
      {
        for (size_t j = i + 1; j < len; j++)
        {
          if ( (points[i].distanceTo(points[j])) <= distance)
          {
            AcDbLine *pLine = new AcDbLine(points[i],points[j]);
            if(pRec->appendAcDbEntity(pLine) != Acad::eOk)
              delete pLine;
            pLine->close();
          }
        }
      }
      pRec->close();
    }
    QueryPerformanceCounter((LARGE_INTEGER*)&end);
    acutPrintf(_T("\n%g Seconds\n"), (double)(end - start) / (double) freq);
  }

arx for 07
« Last Edit: April 02, 2010, 11:19:40 PM by Daniel »

alanjt

  • Needs a day job
  • Posts: 5352
  • Standby for witty remark...
Thank you all :), Your algorithms are all faster than me.
In my PC, (E7500 CPU, 4G Mem)
Ronjonp 's code, about 34.812967s
Lee Mac's code, about 40.733018s
Gile's code 42.111982s
alanjt's code is about 50s
Mine first code is very slow.   :-D, (could you tell me which cause the very slow)
Not very accurate test.

And I hope to see pkohut's super code :)

I am thinking about the "Divide and conquer algorithm", but I need some time to sort the point in some manner, and not so sure about whether it can be improved, because the sort also need some times.

That's interesting that Gile's is 8 seconds faster. Unless I'm overlooking something, we used the same method (comparing with my _connect2), except that I had it broken up into 3 separate subs. I didn't think it would effect performance, but I guess it does. Food for thought.
Civil 3D 2019 ~ Windohz 7 64bit
Dropbox

alanjt

  • Needs a day job
  • Posts: 5352
  • Standby for witty remark...
brute force with arx

doit
GetDist: 70
1.036883 Seconds

arx for 07

That's incredible Daniel! I have got to put Lisp down and learn C.
Civil 3D 2019 ~ Windohz 7 64bit
Dropbox

It's Alive!

  • Retired
  • Needs a day job
  • Posts: 8662
  • AKA Daniel
I'm guessing Paul is pulling apart his triangle routine.. right?

Kerry

  • Mesozoic relic
  • Seagull
  • Posts: 11654
  • class keyThumper<T>:ILazy<T>
brute force with arx

doit
GetDist: 70
1.036883 Seconds

arx for 07

That's incredible Daniel! I have got to put Lisp down and learn C.

:)

Lisp can't be 'put down' it just refuses to die  :-D

I agree about the idea of learning C, C++.  Great if we can make the time and if the wet stuff is still receptive.
kdub, kdub_nz in other timelines.
Perfection is not optional.
Everything will work just as you expect it to, unless your expectations are incorrect.
Discipline: None at all.

pkohut

  • Guest
The 80/20 rule kicked in and it took me a bit longer to write than estimated.

My ARX entry is attached.  ConnectDots2004.arx is for Autocad 2004-2006, and ConnectDots2007 is for 2007-2008(9?).  Just load the appropriate version and at the command line type
Code: [Select]
connectdotsYou will them be prompted for the maximum segment length to connect, and to select the points to operate on.

Here the processing times (counted in my head) -
10,000 points, uh' its done.
100,000 points, little less then 2 seconds.
500,000 points, about 15 seconds.

The heart of the program is the KD-Tree data structure which holds point data selected by the user. Given any point, we can ask the structure to return a set of points within a radius of the source point. Line segments can then be created from each point in the set to the source point. Doing this for all the points in the collection, results in at best 2 line segments for each qualified set of points. There are a number of methods to find duplicate line segments, but they all have some performance cost associated with them. Instead I set a flag for the source point, marking it as already processed and don't have to worry about duplicate segments.

Points of interest in the code (since there is a lot of code)-
Create an instance of CKdTree and initialize it to the number of dimensions. kdTree is set to 3 dimensions, though for the purposes of this program 2 would have been fine since the Z is ignored in all the calculates.
Code: [Select]
CKdTree kdTree;
kdTree.Create(3);

Add a node to the tree. The forth parameter is of type void *, and is for user defined data. In this case I'm using it as a flag to determine if the node has been processed. Initial value is NULL or false.
Code: [Select]
kdTree.Insert(pt3d.x, pt3d.y, pt3d.z, NULL);
Query for all the gathered points (need to add another public function to CKdTree to get the head position).
Get a pointer to the X, Y, Z data
Query for another result set around the X, Y, Z, and radius.
Code: [Select]
RESULT_SET * pResAllPoints = kdTree.NearestRange(0.0, 0.0, 0.0, HUGE_VAL);
...
double * pPos = pResAllPoints->riter->item->pos;
RESULT_SET * pRes = kdTree.NearestRange(pPos[0], pPos[1], 0.0, dRadius);

The KD-Tree code I got off the internet about 5 years ago, though I've not been able to find it again and it had no copyright notice or name associated with it. Any changes I did to the code was so it would compile in C++, IIRC work with doubles instead of floats, and a couple helper functions.  I've used the code in about 7 projects, and it's always fun when another opportunity arises to use it again.

The rest of the code should be straight forward, if not a bit messy.  :-)

pkohut

  • Guest
Have an hour or so of coding left to do.
Man, one really pays for performance. I'm not knocking C/Arx, just taken aback at the amount of time.

Your right, and this is even repurposing code I've already written. The solution though is non-trivial, but screams. Tested with a million points and it's under a second. Still coding though so back to work  :lmao:

When I posted that message I misread the number of points in the selection set. I should have wrote 100,000. sry

It's Alive!

  • Retired
  • Needs a day job
  • Posts: 8662
  • AKA Daniel
blistering fast  8-)

qjchen

  • Bull Frog
  • Posts: 285
  • Best wishes to all
Thank you Daniel for your arx, but I have only 2004 for this PC, I will test your codes later.

And thank you pkohut, Your arx code is very very fast :), I will learn your article later.

After 1 hour programming, I wrote this first lisp revision.

It is about  3.2s for 10000 points, and 6.5s for 20000 points and 17.233984s for 100000 points.

I use many idea of ronjonp, Lee mac ,gile and alanjt, and then, I divide the pointlist to many sub pointlist according to the x coordinate

((points  xmin<x<xmin+70) (points  xmin+70<x<xmin+140) .........)

Then I connect the shortline in each sub pointlist.

Second, I connect the shortline between the adjacent pointlist.

Maybe it can be called "Divide and conquer" method.

In the next revision, I think I should divide the pointlist to x y grid, maybe it will more quickly.


Code: [Select]

(vl-load-com)
;;;start of Time, from theswamp
(defun startTimer ()
  (setq time (getvar "DATE"))
)
;;;end of Time, from theswamp
(defun endTimer (/ seconds totaltime)
  (setq time (- (getvar "DATE") time)
seconds (* 86400.0 (- time (fix time)))
  )
  (setq totaltime (rtos seconds 2 6))
  (princ (strcat "\nTimed:" (rtos seconds 2 6)))
)
;;entmake line
(defun entmakeline (pt1 pt2 layer)
  (entmake (list (cons 0 "LINE") ;***
(cons 6 "BYLAYER") (cons 8 layer) (cons 10 pt1) ;***
(cons 11 pt2) ;***
(cons 39 0.0) (cons 62 256) (cons 210 (list 0.0 0.0 1.0))
   )
  )
)
;;The idea of gile
(defun getpointscor ( / ent n plst ss)
  (if (setq n -1
    ss (ssget '((0 . "POINT")))
      )
    (while (setq ent (ssname ss (setq n (1+ n))))
      (setq plst (cons (cdr (assoc 10 (entget ent))) plst))
    )
  )
)
(defun dividex (plst dis / a b div end index len lst1 plstmp res start xend
     xstart
       )
  (setq start (car plst)
xstart (car start)
end (last plst)
len (- (car end) (car start))
div (+ (fix (/ len dis)) 1)
plstmp plst
res nil
  )
  (repeat div 
    (grdraw (list xstart 0 0) (list xstart 10000 0) 8)
    (setq xend (+ xstart dis))
    (setq index T
  lst1 nil
    )
    (while (and
     plstmp
     index
   )
      (if (< (setq a (car plstmp)
   b (car a)
     )
     xend
  )
(setq lst1 (append
     lst1
     (list a)
   )
      plstmp (cdr plstmp)
)
(setq index nil)
      )
    )
    (setq res (append
res
(list lst1)
      )
    )
    (setq xstart (+ xstart dis))
  )
  res
)
;;;draw short line in each x region
;;;The idea of Leemac and gile
(defun drawshortline (plst limdis / a b tmp)
  (while (setq a (car plst))
    (setq plst (cdr plst)
  tmp plst
    )
    (while (setq b (car tmp))
      (setq tmp (cdr tmp))
      (if (< (distance a b) limdis)
(entmakeline a b "temp")
      )
    )
  )
)
;;;draw short line in each two Neighbor x region
(defun drawshortline12 (plst1 plst2 limdis / a b x)
  (while (setq a (car plst1))
    (foreach x plst2
      (if (< (distance a x) limdis)
(entmakeline a x "temp")
      )
    )
    (setq plst1 (cdr plst1))
  )
)


;;Connect short line between each two points
;;by qjchen@gmail.com
;; Thanks to ronjonp, Lee Mac, and Gile
;;Version 1.0
(defun c:gap ( / a b limdis plst plst1)  ;(mylayer "temp")
  (starttimer)
  (setq limdis 70)
  ;;The idea of gile
  (setq plst (getpointscor))
  ;;The idea of ronjonp
  (setq plst (vl-sort plst (function (lambda (x1 x2)(< (car x1) (car x2))))))
  ;;Divde points list to several region
  ;;such like ((0, contain x=0~70),(1, contain x=70~140)........)
  (setq plst1 (dividex plst limdis))
  ;;Draw shortline in each x region ,such like (0, contain x=0~70))
  (foreach x plst1
    (drawshortline x limdis)
  )
  ;;Draw shortline between two Neighbor x region, such between (0, contain x=0~70),(1, contain x=70~140))
  (while (setq a (car plst1) b (cadr plst1) plst1 (cdr plst1))
    (drawshortline12 a b limdis)
  )
  (endtimer)
)


« Last Edit: April 03, 2010, 11:00:54 AM by yuanqiu »
http://qjchen.mjtd.com
My blog http://chenqj.blogspot.com (Chinese, can be translate into English)

alanjt

  • Needs a day job
  • Posts: 5352
  • Standby for witty remark...
brute force with arx

doit
GetDist: 70
1.036883 Seconds

arx for 07

That's incredible Daniel! I have got to put Lisp down and learn C.

:)

Lisp can't be 'put down' it just refuses to die  :-D

I agree about the idea of learning C, C++.  Great if we can make the time and if the wet stuff is still receptive.

I just meant to put it down long enough to start with a new language. It's far too useful to abandon.
Civil 3D 2019 ~ Windohz 7 64bit
Dropbox

ElpanovEvgeniy

  • Water Moccasin
  • Posts: 1569
  • Moscow (Russia)
my variant:

Code: [Select]
(defun el (d / A L LL TI X)
 ;;(el 70.)
 (setq ti (car (_VL-TIMES))
       a  -1
       ll (ssget "_x" '((0 . "point")))
 )
 (while (setq x (ssname ll (setq a (1+ a))))
  (setq l (cons (cdr (assoc 10 (entget x))) l))
 )
 (setq l (vl-sort l (function (lambda (a b) (<= (car a) (car b))))))
 (while (setq a (car l)
              l (cdr l)
        )
  (setq x  (car a)
        ll l
  )
  (while (and ll (<= (- (caar ll) x) d))
   (if (and (equal (cdar ll) (cdr a) d) (<= (distance a (car ll)) d))
    (entmakex (list '(0 . "LINE") (cons 10 a) (cons 11 (car ll))))
   )
   (setq ll (cdr ll))
  )
 )
 (princ (strcat "\n " (rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4) "sec."))
 (princ)
)

It's Alive!

  • Retired
  • Needs a day job
  • Posts: 8662
  • AKA Daniel
Command: doit
497.0208 Milliseconds 8-)

Code: [Select]
using System;
using System.Collections.Generic;
using System.Text;
using System.Threading;
using System.Threading.Tasks;
using System.Diagnostics;

using Autodesk.AutoCAD.GraphicsInterface;
using Autodesk.AutoCAD.Runtime;
using Autodesk.AutoCAD.Windows;
using Autodesk.AutoCAD.Geometry;
using Autodesk.AutoCAD.EditorInput;
using Autodesk.AutoCAD.LayerManager;
using Autodesk.AutoCAD.ApplicationServices;
using Autodesk.AutoCAD.DatabaseServices;
using Autodesk.AutoCAD.Colors;
using System.Threading.Tasks;


using AcAp = Autodesk.AutoCAD.ApplicationServices;
using AcEd = Autodesk.AutoCAD.EditorInput;
using AcGe = Autodesk.AutoCAD.Geometry;
using AcRx = Autodesk.AutoCAD.Runtime;
using AcDb = Autodesk.AutoCAD.DatabaseServices;
using AcWd = Autodesk.AutoCAD.Windows;

namespace Test
{
  public class ExecMethod
  {
    [CommandMethod("doit")]
    public static void MyCommand()
    {

      Stopwatch stopWatch = new Stopwatch();
      stopWatch.Start();

      SynchronizedCollection<Line>
        lines = new SynchronizedCollection<Line>();

      List<Point3d> points = getPoints();
      int cnt = points.Count;


      Parallel.For(0, cnt, delegate(int i)
      {
        Point3d pnt = points[i];
        for (int j = i + 1; j < cnt; j++)
        {
          if (pnt.DistanceTo(points[j]) <= 70.0)
          {
            lines.Add(new Line(pnt, points[j]));
          }
        }
      });
      addtoDb(lines);

      stopWatch.Stop();
      TimeSpan ts = stopWatch.Elapsed;
      AcAp.Application.DocumentManager.MdiActiveDocument.
        Editor.WriteMessage("\n{0} Seconds", ts.TotalMilliseconds);
    }

    static void addtoDb(SynchronizedCollection<Line> lines)
    {
      List<Point3d> points = new List<Point3d>();
      Database database = HostApplicationServices.WorkingDatabase;
      using (Transaction transaction = database.TransactionManager.StartTransaction())
      {
        BlockTableRecord record =
          transaction.GetObject(
            database.CurrentSpaceId, OpenMode.ForWrite) as BlockTableRecord;

        foreach (Line line in lines)
        {
          record.AppendEntity(line);
          transaction.AddNewlyCreatedDBObject(line, true);
        }
        transaction.Commit();
      }
    }


    static List<Point3d> getPoints()
    {
      List<Point3d> points = new List<Point3d>();
      Database database = HostApplicationServices.WorkingDatabase;
      using (Transaction transaction = database.TransactionManager.StartTransaction())
      {
        BlockTableRecord record =
          transaction.GetObject(
            database.CurrentSpaceId, OpenMode.ForRead) as BlockTableRecord;

        foreach (ObjectId id in record)
        {
          DBPoint point = transaction.GetObject(id, OpenMode.ForRead) as DBPoint;
          if (point != null)
          {
            points.Add(point.Position);
          }
        }
        transaction.Commit();
      }
      return points;
    }
  }
}

It's Alive!

  • Retired
  • Needs a day job
  • Posts: 8662
  • AKA Daniel

alanjt

  • Needs a day job
  • Posts: 5352
  • Standby for witty remark...
Civil 3D 2019 ~ Windohz 7 64bit
Dropbox

It's Alive!

  • Retired
  • Needs a day job
  • Posts: 8662
  • AKA Daniel
my variant:

Code: [Select]
(defun el (d / A L LL TI X)
 ;;(el 70.)
 (setq ti (car (_VL-TIMES))
       a  -1
       ll (ssget "_x" '((0 . "point")))
 )
 (while (setq x (ssname ll (setq a (1+ a))))
  (setq l (cons (cdr (assoc 10 (entget x))) l))
 )
 (setq l (vl-sort l (function (lambda (a b) (<= (car a) (car b))))))
 (while (setq a (car l)
              l (cdr l)
        )
  (setq x  (car a)
        ll l
  )
  (while (and ll (<= (- (caar ll) x) d))
   (if (and (equal (cdar ll) (cdr a) d) (<= (distance a (car ll)) d))
    (entmakex (list '(0 . "LINE") (cons 10 a) (cons 11 (car ll))))
   )
   (setq ll (cdr ll))
  )
 )
 (princ (strcat "\n " (rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4) "sec."))
 (princ)
)

you're good  8-)