Calculate outer paths, given centre path and point radius

Just starting out? Need help? Post your questions and find answers here.
Seymour Clufley
Addict
Addict
Posts: 1233
Joined: Wed Feb 28, 2007 9:13 am
Location: London

Calculate outer paths, given centre path and point radius

Post by Seymour Clufley »

My code takes a set of points (that form a zig-zaggy line) and attempts to find the outer line, on either side of that centre line, given a width - and the width can vary at each point.

First, you give it the points for the centre line:
Image

Then the radius for each point:
Image

And it calculates the "outer lines", returning the equivalent point set for each one:
Image

It works well most of the time, but then something like this happens:
Image

I have written code to eliminate these errors after the procedure - it loops through the lines, checking whether lines N and N+2 intercept each other - and it works, but I think there is probably a way to eliminate these errors before they are made. If anyone can help with this, I'd be grateful.

Here is the code:

Code: Select all


Structure PointD
  x.d
  y.d
EndStructure
Structure LineD
	pnt1.PointD
	pnt2.PointD
	deg.d
	length.d
EndStructure



Macro StandardReportingWindowStuff(win)
   escapekey = 1
   returnkey = 2
   AddKeyboardShortcut(win,#PB_Shortcut_Escape,escapekey)
   AddKeyboardShortcut(win,#PB_Shortcut_Return,returnkey)
   
   MoveWindowToTop(win,#True)
   
   
   Repeat
      we = WindowEvent()
      If we
         If we=#PB_Event_Menu
            Break
         EndIf
      Else
         Delay(10)
      EndIf
   ForEver
   
   CloseWindow(win)
   
EndMacro

Procedure.b MoveWindowToTop(win.i,stayontop.b)
   If Not WindowID(win) : MessageRequester("Error","CAN'T FIND WINDOW ID") : ProcedureReturn #False : EndIf
   If stayontop
      SetWindowPos_(WindowID(win),#HWND_TOPMOST,0,0,0,0,#SWP_NOMOVE|#SWP_NOSIZE)    ; gets it on top to stay 
   Else
      SetWindowPos_(WindowID(win),#HWND_NOTOPMOST,0,0,0,0,#SWP_NOMOVE|#SWP_NOSIZE)  ; now it is on top, allow it to go behind
   EndIf
   ProcedureReturn #True
EndProcedure

Procedure.b ResizeToFitInsideFrame(iw.d,ih.d,fw.d,fh.d,*ns.PointD)
   
   *ns\x = iw
   *ns\y = ih
   
   While *ns\x>fw Or *ns\y>fh
     *ns\x * 0.99
     *ns\y * 0.99
   Wend
   ProcedureReturn #True
   
EndProcedure

Procedure RI(img.i,title.s="",bgclr.i=#Green,fit_on_screen.b=#True,timelimit.i=0,free.b=#False)
  
  If Not IsImage(img)
    MessageRequester("Error","PROC RI (ReportImage):"+c13+"img is not an image!")
    ProcedureReturn
  EndIf
  
  iw.d = ImageWidth(img)
  ih.d = ImageHeight(img)
  
  simg = CreateImage(#PB_Any,iw,ih,32,bgclr)
  StartDrawing(ImageOutput(simg))
  DrawAlphaImage(ImageID(img),0,0)
  StopDrawing()
  
  
  If title=""
    title = "Report Image"
  EndIf
  
  win = OpenWindow(#PB_Any,0,0,iw,ih,title,#PB_Window_BorderLess|#PB_Window_ScreenCentered)
  
  imgad = ImageGadget(#PB_Any,0,0,iw,ih,ImageID(simg))
  
  StandardReportingWindowStuff(win)
  
  FreeImage(simg)
  
EndProcedure



Procedure.d Difference(a.d,b.d)
  
  If a=b
      ProcedureReturn 0
  EndIf
  If a>b
  	ProcedureReturn a-b
  EndIf
  ProcedureReturn b-a
  
EndProcedure

Procedure.d DistanceBetweenTwoPoints(*a.PointD,*b.PointD)
	xdif.d = Difference(*a\x,*b\x)
	ydif.d = Difference(*a\y,*b\y)
	ProcedureReturn Sqr((xdif*xdif)+(ydif*ydif))
EndProcedure

Procedure.d DegreeAngleBetweenTwoPoints(*o.PointD,*b.PointD)
  calcAngle.d = Degree(ATan2(*b\x-*o\x,*b\y-*o\y))
  If calcAngle<0
    calcAngle = 360-Abs(calcAngle)
  EndIf
  ProcedureReturn calcAngle
EndProcedure


Procedure.b RadianCoordsFromPoint(*base.PointD,radia.d,distance.d,*c.PointD)
	*c\x = (Cos(radia)*distance)+*base\x
	*c\y = (Sin(radia)*distance)+*base\y
EndProcedure

Procedure.b DegreeCoordsFromPoint(*base.PointD,degrees.d,distance.d,*c.PointD)
	RadianCoordsFromPoint(*base,Radian(degrees),distance,*c)
EndProcedure


Procedure.b PointsVirtuallyIdentical(*p1.PointD,*p2.PointD)
	ProcedureReturn Bool( Difference(*p1\x,*p2\x)<1 And Difference(*p1\y,*p2\y)<1 )
EndProcedure


Procedure.b TwoInfiniteLinesIntersection(*l1.LineD,*l2.LineD,*i.PointD)
	
	If PointsVirtuallyIdentical(*l1\pnt2,*l2\pnt1)
		CopyStructure(*l1\pnt2,*i,PointD)
		ProcedureReturn #True
	EndIf
	
	a1.d = *l1\pnt2\y-*l1\pnt1\y
	b1.d = *l1\pnt1\x-*l1\pnt2\x
	c1.d = ( *l1\pnt2\x * *l1\pnt1\y ) - ( *l1\pnt1\x * *l1\pnt2\y )
	
	a2.d = *l2\pnt2\y-*l2\pnt1\y
	b2.d = *l2\pnt1\x-*l2\pnt2\x
	c2.d = ( *l2\pnt2\x * *l2\pnt1\y ) - ( *l2\pnt1\x * *l2\pnt2\y )
	
	denom.d = (a1*b2) - (a2*b1)
	*i\x = (b1*c2 - b2*c1)/denom
	*i\y = (a2*c1 - a1*c2)/denom
	ProcedureReturn #True
	
EndProcedure


Procedure.b OuterLines(pnts.i,Array pnt.PointD(1),Array rad.d(1),Array outerside1.PointD(1),Array outerside2.PointD(1))
  
  Dim leftside.LineD(pnts)
  Dim rightside.LineD(pnts)
  
  For p = 1 To pnts-1
    deg.d = DegreeAngleBetweenTwoPoints(@pnt(p),@pnt(p+1))
    dist.d = DistanceBetweenTwoPoints(@pnt(p),@pnt(p+1))
    
    DegreeCoordsFromPoint(@pnt(p),deg-90,rad(p),@leftside(p)\pnt1)
    ;Circle(leftside(p)\pnt1\x,leftside(p)\pnt1\y,4,#Yellow)
    DegreeCoordsFromPoint(@pnt(p+1),deg-90,rad(p+1),@leftside(p)\pnt2)
    ;Circle(leftside(p)\pnt2\x,leftside(p)\pnt2\y,4,#Blue)
    
    DegreeCoordsFromPoint(@pnt(p),deg+90,rad(p),@rightside(p)\pnt1)
    DegreeCoordsFromPoint(@pnt(p+1),deg+90,rad(p+1),@rightside(p)\pnt2)
    ;Circle(rightside1\x,rightside1\y,4,#Red)
    
  Next p
  
  Dim outerside1(pnts)
  Dim outerside2(pnts)
  ; first point
  CopyStructure(@leftside(1)\pnt1,@outerside1(1),PointD)
  CopyStructure(@rightside(1)\pnt1,@outerside2(1),PointD)
  ; the others
  For p = 1 To pnts-1
    TwoInfiniteLinesIntersection(@leftside(p),@leftside(p+1),@outerside1(p+1))
    TwoInfiniteLinesIntersection(@rightside(p),@rightside(p+1),@outerside2(p+1))
  Next p
  ; the final point
  deg.d = DegreeAngleBetweenTwoPoints(@pnt(pnts-1),@pnt(pnts))
  dist.d = DistanceBetweenTwoPoints(@pnt(pnts-1),@pnt(pnts))
  DegreeCoordsFromPoint(@pnt(pnts),deg-90,rad(pnts),@outerside1(p))
  DegreeCoordsFromPoint(@pnt(pnts),deg+90,rad(pnts),@outerside2(p))
  
EndProcedure




iw.d = 2000
ih.d = 1100

Dim pnt.PointD(30)
Dim rad.d(30)
dx.d = 100
p=0
  While dx<iw-100
    p+1
    pnt(p)\x = dx
		pnt(p)\y = (ih/2)+Random(ih/3)-Random(ih/3)
		rad(p) = 10+Random(47)
  	dx+50+Random(150)
Wend
pnts = p


; SET 1: "triangle" error
; pnts = 4
; pnt(1)\x=150 : pnt(1)\y=50
; rad(1)=30
; pnt(2)\x=255 : pnt(2)\y=115
; rad(2)=35
; pnt(3)\x=250 : pnt(3)\y=130
; rad(3)=20
; pnt(4)\x=150 : pnt(4)\y=380


; SET 2: another "triangle" error
; pnts = 6
; pnt(1)\x=750 : pnt(1)\y=50
; rad(1)=33
; pnt(2)\x=767.9718826698 : pnt(2)\y=238
; rad(2)=28
; pnt(3)\x=701.6245164548 : pnt(3)\y=282
; rad(3)=24
; pnt(4)\x=718.4686470836 : pnt(4)\y=290
; rad(4)=19
; pnt(5)\x=769.8000570659 : pnt(5)\y=326
; rad(5)=14
; pnt(6)\x=750 : pnt(6)\y=399
; rad(6)=9




img = CreateImage(#PB_Any,iw,ih,32)
StartDrawing(ImageOutput(img))
DrawingMode(#PB_2DDrawing_AlphaBlend)

; show the points...
For p = 1 To pnts
  Circle(pnt(p)\x,pnt(p)\y,rad(p),RGBA(255,255,0,48))
Next p

; draw the centre line...
For p = 1 To pnts-1
  LineXY(pnt(p)\x,pnt(p)\y,pnt(p+1)\x,pnt(p+1)\y,RGBA(255,255,255,160))
Next p

Dim outer1.PointD(pnts)
Dim outer2.PointD(pnts)
OuterLines(pnts,pnt(),rad(),outer1(),outer2())
For p = 1 To pnts-1
  ;Circle(outer1(p)\x,outer1(p)\y,1.5,RGBA(255,255,255,255))
  ;Circle(outer2(p)\x,outer2(p)\y,1.5,RGBA(255,255,0,255))
  LineXY(outer1(p)\x,outer1(p)\y,outer1(p+1)\x,outer1(p+1)\y,RGBA(0,255,0,255))
  LineXY(outer2(p)\x,outer2(p)\y,outer2(p+1)\x,outer2(p+1)\y,RGBA(0,255,0,255))
Next p

StopDrawing()
RI(img)
FreeImage(img)
It will generate a random path each time it is run. Two "problematic" paths are included, commented out.
JACK WEBB: "Coding in C is like sculpting a statue using only sandpaper. You can do it, but the result wouldn't be any better. So why bother? Just use the right tools and get the job done."
IdeasVacuum
Always Here
Always Here
Posts: 6425
Joined: Fri Oct 23, 2009 2:33 am
Location: Wales, UK
Contact:

Re: Calculate outer paths, given centre path and point radiu

Post by IdeasVacuum »

Hi Seymour

This is a classic "bow tie" issue in CAD. I think your solution is already good - just don't draw the outer path lines until after the clean-up routine has performed it's magic.

Would the Vector Drawing Lib make life easier?
IdeasVacuum
If it sounds simple, you have not grasped the complexity.
Post Reply