FB II Compiler

PG PRO

Debugging

Memory

System

Mathematics

Resources

Disk I/O

Windows

Controls

Menus

Mouse

Keyboard

Text

Fonts

Drawing

Sound

Clipboard

Printing

Communication

ASM

Made with FB

MATHEMATICS

Solve a geometric problem


I wonder if anyone anyone can help me with this little dilemma that I must solve to get an image processing program to work.
My program finds 3 points in a picture which I have the x and y coordinates of: Ax,Ay, Bx,By and Cx,Cy. I want to figure out where a fourth point would be. I want to find the point on a line that passes through A and B that is closest to C. If anyone knows the method to do this the help would be greatly appreciated. I knew I should have paid more attention in algebra class in high school.

Joe Lertola


FBII I'm learning! Analytic geometry I think I know.
The point on a line closest to a given point is found by dropping a perpendicular from the point to the line.
While you could develop an explicit expression for DX and DY, it's probably best to do it in pieces.
The slope of the line connecting A & B, call it M, is (AY-BY)/(AX-BX) so...
The equation of the line connecting A & B is y-AY M*(x-AX) or
y = (M*x) + (AY-M*AX)
The slope of the line perpendicular to the line connecting A & B is -1/M, call it N, so...
The equation of the line through C perpendicular to A&B is y-CY = N*(x-CX) or
y = (N*x) + (CY-N*CX)
So you have two linear equations in two unknowns. If you equate AY-M*AX to P and CY-N*CX to Q then the two equations are:
y = Mx + P
y = Nx + Q
The solution to this system of equations, which will be DX and DY, is:
DX = (Q-P)/(M-N)
DY = (QM-NP)/(M-N)
You have to check for special cases like horizontal and vertical lines.
I hope this helps. I hope I didn't make an algebra error.

Jeff Schwartz


The equation of the line containing the points A and B is y = mx + b, where m = ((y(B)-y(A))/((x(B)-x(A)) and b = y(A) - mx(A)
The distance you want to minimize is the distance from point C to a point on the line. That distance (or the square of it) is
d = (x(C)-x)^2 + (y(C)-y)^2
Taking the derivative and setting it to zero and solving for x that gives the minimum distance gives
x = (x(C)-m(b-y(C)))/(1+m^2)
y is then gotten from y = mx +b
x and y above are the coordinates of the point you want assuming I made no errors. It should always be valid for any three points that don't fall on a line.

Eric


OK, here goes...
You want the coordinates of D (Dx, Dy) which is closest to the line formed by A and B (let's call it AB).
This point will be on AB and on the line perpendicular to AB through C (the perpendicular distance is always the shortest).
Let M = (By-Ay)/(Bx-Ax). Then M is the slope of AB and -1/M is the slope of CD. Finding the coordinates of D goes like this...
Dx = ((Cx+((M^2)*Ax)-(M*Ay)+(M*Cy))/((M^2)+1) and Dy = (M*(Dx-Ax))+Ay
While there's a lot of algebra involved, the subject you missed is analytic geometry.

Charlie Dickman


My thanks to Jeff Schwartz, Eric and Charlie Dickman for your help with this problem. I can't say that I follow any of the explanations but as long as I have some code that works I am happy.

Jeff, This is how I interpret your equations as FutureBasic code. It seems to work perfectly even though I can't quite seem to follow your explanation.

WINDOW 1
Ax=0:Ay=0
Bx=5:By=5
Cx=3:Cy=1
M#=(Ay-By)/(Ax-Bx)
N#=-1/M#
P#=Ay-M#*Ax
Q#=Cy-N#*Cx
Dx= (Q#-P#)/(M#-N#)
Dy= (Q#*M#-N#*P#)/(M#-N#)
PRINT Dx,Dy

Eric, This is how I interpret your equations. It does not seem to produce correct results. I may be misinterpreting what you wrote.

m#=(By-Ay)/(Bx-Ax)
b# = Ay - m# * Ax
Dx=(Cx-m#*(b#-Cy))/1+(m#^2)
Dy=m#*Dx+b#
PRINT Dx,Dy

Charles, This is your code, and it works great. It is about 30% faster than Jeff's. Do you see any way to speed it up even more? I intend to use this as part of a program that process big image files and this will be applied to every pixel. I need to optimize it for speed.

M# = (By-Ay)/(Bx-Ax)
Dx = (Cx+((M#^2)*Ax)-(M#*Ay)+(M#*Cy))/((M#^2)+1)
Dy = (M#*(Dx-Ax))+Ay
PRINT Dx,Dy

Joe Lertola


You have to be careful to test before computing M# that Bx<>Ax in this case the solution is evident (Dx=Ax) but the general algorithm fails.
To speed up the computation, I suggest

long if Ax<>BX
  M# = (By-Ay)/(Bx-Ax)
  M2#= M#*M#
  Dx = (Cx+(M2#*Ax)-(M#*Ay)+(M#*Cy))/(M2#+1)
  Dy = (M#*(Dx-Ax))+Ay
xelse
  Dx=Ax
  Dy=Cy
end if

Guy Bonnot


The equations needed a FB tweaking since parentheses and other operations convert values to integers in the formulas.
I have provided some code that plots the solution. There are 2 FNs. The first is for use with integer values (not used in my code)
The second is for use with floating values (used in the code).

CLEAR LOCAL 'Uses integer points
  DIM M#
  DIM M2#
  DIM Dx#
  DIM Dy#
LOCAL FN FindNearPtinLine(Ax,Ay,Bx,By,Cx,Cy)
  LONG IF Ax<>Bx
    M# = (By-Ay)/(Bx-Ax)
    M2#= M#*M#
    Dx# = (Cx+M2#*Ax-M#*Ay+M#*Cy)\(M2#+1.0)
    Dy# = M#*(Dx#-Ax)+Ay
  XELSE
    Dx#=Ax
    Dy#=Cy
  END IF
  COLOR =_zRed
  CIRCLE FILL Dx#,Dy#,2
  COLOR =_zBlack
END FN

CLEAR LOCAL 'Uses float points
  DIM M#
  DIM M2#
  DIM Dx#
  DIM Dy#
LOCAL FN FindNearFloatPtinLine(Ax#,Ay#,Bx#,By#,Cx#,Cy#)
  LONG IF Ax#<>Bx#
    M# = (By#-Ay#)/(Bx#-Ax#)
    M2#= M#*M#
    Dx# = (Cx#+(M2#*Ax#)-(M#*Ay#)+(M#*Cy#))/(M2#+1.0)
    Dy# = (M#*(Dx#-Ax#))+Ay#
  XELSE
    Dx#=Ax#
    Dy#=Cy#
  END IF
  COLOR =_zRed
  CIRCLE FILL Dx#,Dy#,3
  PLOT Cx#,Cy# TO Dx#,Dy#
  COLOR =_zBlack
END FN

WINDOW 1
Ax#=5
Ay#=20
Bx#=100
By#=120
Cx#=50
Cy#=100
FOR y= 1 TO 5
By#=By#-20
CIRCLE Ax#,Ay#,3
CIRCLE Bx#,By#,3
PLOT Ax#,Ay# TO Bx#,By#
COLOR =_zRed
CIRCLE FILL Cx#,Cy#,3
FN FindNearFloatPtinLine(Ax#,Ay#,Bx#,By#,Cx#,Cy#)
NEXT y
DO
UNTIL FN BUTTON

Pierre A. Zippi


Yes, there are ways to improve the code...
calculate M#^2 only once and save it in a variable, say M2#, and then use it in the expression for Dx in place of (M#2).
Where you have ...-(M#*Ay)+(M#*Cy)... factor out the M# like this... M#*(Cy-Ay) so the expression for Dx is Dx = (Cx+(M2*Ax)+(M#*(Cy-Ay)))/(M2#+1). This reduces the number of multiplications by 3.
Also, you must special case the situation when AB is a vertical line (Bx-Ax = 0) or you will divide by 0 when computing M and this is a no-no. In this case Dx = Bx (or Dx = Ax) and Dy = Cy. You can also special case the case of AB being a horizontal line (By - Ay = 0). This case reduces to Dx = Cx and Dy = Ay (or Dy = By).

Charlie Dickman


Speaking of streamlining equations for speed...
What is the fastest method for rotating an x-y plot.
Lets say that I have a standard horizontal-vertical x-y plot and and I now wish to re-plot the same relative x-y plot but with the x and y axes rotated at some angle.
When trig functions get involved things always slow down.
To compound matters what I REALLY want to do is adjust the position and angle of the x and y axes (always perpendicular) with every point plotted. So there is no guarantee that the rotation angle will remain constant from point to point.
Am I stuck with rotating every point with trig functions?
If so, what is the speediest equation?

Pierre A. Zippi


I have had some experience speeding things up. One of the tricks I use is to never use floating point numbers inside of loops. The recent example you posted that finds a point on a line took about 47 ticks to run 10000 times. I was able to get this down to 13 ticks by removing all the floating point calculations. For example instead of:

M# = (By-Ay)/(Bx-Ax)

I did:

Mtop&=(By-Ay)
Mbot&=(Bx-Ax)

And then any where in the calculation where I should have multiplyed by M#, I multipiled by Mtop& then divided by Mbot&. Instead of:

Dy = (M#*(Dx-Ax))+Ay

I did:

Dy = (Mtop&*(Dx-Ax)/Mbot&)+Ay

(By the way Pierre, thanks for the example you posted. It was a big help. Also thanks to Guy Bonnot, Charlie Dickman, Jeff Schwartz and Eric)
I had a program that I wrote several years ago that cropped large geographic data bases. The first version used floating point numbers in its calculations.
A typical data base took about 4 hours to be processed by this program. When I optimized the program for speed by removing all the floating point calculations the program took only 2 minutes to crop a typical data base.
Another technique to use is to precalculate your function values and store them in a array. For example build a table for sin values like this (not tested):

dim sinVal&(360)
for x =0 to 360
  Pi#=3.1415926535898
  radian#=((2*Pi#)/360)*x
  sinFloat#=sin(radian#)
  sinVal&(x)=sinFloat#*1000 '<- multiply by some appropriate scaling number so you can store you number as a long integer
next x

Then, when you are in your main loop, look up the sin of an angle from the array, which is much faster than calculating it. Divide the scaling number out of your answer when you are done with your calculation.
If you use just integer math and do not use math functions in your main loop, your program will scream along. It will seem as if you have a PPC native app.

Joe Lertola