Author Topic: Ellipse Intersect Line  (Read 520 times)

0 Members and 1 Guest are viewing this topic.

Offline STxAxTIC

  • Library Staff
  • Forum Resident
  • Posts: 1062
  • TOXIC
Ellipse Intersect Line
« on: February 01, 2020, 04:22:30 PM »
In response to

https://www.qb64.org/forum/index.php?topic=2132

... all I could think is "why stop at circles when you can do ellipses?"

So this program calculates and plots the intersection points between any ellipse and any line. It's not polished whatsoever - you may break something if you try and tweak it. I didn't spend the time to make this more software-like. That can be on you guys if needed. Here's the hard part though:

Code: QB64: [Select]
  1. PAINT (1, 1), 15
  2.  
  3. xorig = 0
  4. yorig = 50
  5.  
  6. CALL cline(xorig, yorig, xorig + _WIDTH, yorig, 8)
  7. CALL cline(xorig, yorig, xorig + -_WIDTH, yorig, 8)
  8. CALL cline(xorig, yorig, xorig, yorig + _HEIGHT, 8)
  9. CALL cline(xorig, yorig, xorig, yorig - _HEIGHT, 8)
  10.  
  11. xzoom = 20
  12. yzoom = 20
  13.  
  14. ' Line
  15. b = -2
  16. vy = 1
  17. vx = -.25
  18. v = SQR(vx ^ 2 + vy ^ 2)
  19. vx = vx / v
  20. vy = vy / v
  21. m = vy / vx
  22. FOR alpha = -20 TO 20 STEP .001
  23.     x = alpha * vx
  24.     y = alpha * vy + b
  25.     CALL ccircle(xorig + x * xzoom, yorig + y * yzoom, 1, 1)
  26.  
  27. ' Ellipse
  28. x0 = 2
  29. y0 = -2
  30. ax = 10
  31. ay = 10
  32. bx = -5
  33. by = 5
  34. FOR t = 0 TO 6.284 STEP .001
  35.     x = x0 + ax * COS(t) + bx * SIN(t)
  36.     y = y0 + ay * COS(t) + by * SIN(t)
  37.     CALL ccircle(xorig + x * xzoom, yorig + y * yzoom, 1, 4)
  38.  
  39. ' Intersections
  40. a2 = ax ^ 2 + ay ^ 2
  41. b2 = bx ^ 2 + by ^ 2
  42. av = ax * vx + ay * vy
  43. bv = bx * vx + by * vy
  44. rbx = 0 - x0
  45. rby = b - y0
  46. adbr = ax * rbx + ay * rby
  47. bdbr = bx * rbx + by * rby
  48. aa = av ^ 2 / a2 ^ 2 + bv ^ 2 / b2 ^ 2
  49. bb = 2 * (av * adbr / a2 ^ 2 + bv * bdbr / b2 ^ 2)
  50. cc = adbr ^ 2 / a2 ^ 2 + bdbr ^ 2 / b2 ^ 2 - 1
  51. alpha1 = (-bb + SQR(bb ^ 2 - 4 * aa * cc)) / (2 * aa)
  52. alpha2 = (-bb - SQR(bb ^ 2 - 4 * aa * cc)) / (2 * aa)
  53. x1 = alpha1 * vx
  54. x2 = alpha2 * vx
  55. y1 = alpha1 * vy + b
  56. y2 = alpha2 * vy + b
  57. CALL ccircle(xorig + x1 * xzoom, yorig + y1 * yzoom, 10, 1)
  58. CALL ccircle(xorig + x2 * xzoom, yorig + y2 * yzoom, 10, 1)
  59.  
  60.  
  61. SUB cline (x1, y1, x2, y2, col)
  62.     LINE (_WIDTH / 2 + x1, -y1 + _HEIGHT / 2)-(_WIDTH / 2 + x2, -y2 + _HEIGHT / 2), col
  63.  
  64. SUB ccircle (x1, y1, r, col)
  65.     CIRCLE (_WIDTH / 2 + x1, -y1 + _HEIGHT / 2), r, col
TOXIC

Marked as best answer by STxAxTIC on February 01, 2020, 08:35:47 PM

Offline STxAxTIC

  • Library Staff
  • Forum Resident
  • Posts: 1062
  • TOXIC
Re: Ellipse Intersect Line
« Reply #1 on: February 01, 2020, 06:33:01 PM »
Alright, I made it into "software"... Now you can change everything on the fly. Let me know if any bugs show up. (Note the line is assumed to have infinite length, it's not a ray.)

Code: QB64: [Select]
  1.  
  2. xorig = 0
  3. yorig = 0
  4.  
  5. CALL cline(xorig, yorig, xorig + _WIDTH, yorig, 8)
  6. CALL cline(xorig, yorig, xorig + -_WIDTH, yorig, 8)
  7. CALL cline(xorig, yorig, xorig, yorig + _HEIGHT, 8)
  8. CALL cline(xorig, yorig, xorig, yorig - _HEIGHT, 8)
  9.  
  10. xzoom = 20
  11. yzoom = 20
  12.  
  13. ' Initialize line
  14. b = -2
  15. d = 2
  16. lineang = .1
  17. vx = COS(lineang)
  18. vy = SIN(lineang)
  19. m = vy / vx
  20.  
  21. ' Initialize ellipse
  22. x0 = 2
  23. y0 = -2
  24. ellipsearg = .2
  25. amag = 10
  26. ax = amag * COS(ellipsearg)
  27. ay = amag * SIN(ellipsearg)
  28. bmag = 5
  29. bx = bmag * COS(ellipsearg + 3.14 / 2)
  30. by = bmag * SIN(ellipsearg + 3.14 / 2)
  31.  
  32.  
  33.         x = _MOUSEX
  34.         y = _MOUSEY
  35.         IF ((x > 0) AND (x < _WIDTH) AND (y > 0) AND (y < _HEIGHT)) THEN
  36.             IF _MOUSEBUTTON(1) THEN
  37.                 x = _MOUSEX
  38.                 y = _MOUSEY
  39.                 x0 = (x - _WIDTH / 2) / xzoom
  40.                 y0 = (-y + _HEIGHT / 2) / yzoom
  41.             END IF
  42.             IF _MOUSEBUTTON(2) THEN
  43.                 x = _MOUSEX
  44.                 y = _MOUSEY
  45.                 d = (x - _WIDTH / 2) / xzoom
  46.                 b = (-y + _HEIGHT / 2) / yzoom
  47.             END IF
  48.             IF _MOUSEWHEEL > 0 THEN
  49.                 lineang = lineang + .01
  50.                 vx = COS(lineang)
  51.                 vy = SIN(lineang)
  52.                 m = vy / vx
  53.             END IF
  54.             IF _MOUSEWHEEL < 0 THEN
  55.                 lineang = lineang - .01
  56.                 vx = COS(lineang)
  57.                 vy = SIN(lineang)
  58.                 m = vy / vx
  59.             END IF
  60.         END IF
  61.     LOOP
  62.  
  63.         CASE 18432
  64.             bmag = bmag + .1
  65.             bx = bmag * COS(ellipsearg + 3.14 / 2)
  66.             by = bmag * SIN(ellipsearg + 3.14 / 2)
  67.         CASE 20480
  68.             bmag = bmag - .1
  69.             bx = bmag * COS(ellipsearg + 3.14 / 2)
  70.             by = bmag * SIN(ellipsearg + 3.14 / 2)
  71.         CASE 19200
  72.             ellipsearg = ellipsearg - .1
  73.             ax = amag * COS(ellipsearg)
  74.             ay = amag * SIN(ellipsearg)
  75.             bx = bmag * COS(ellipsearg + 3.14 / 2)
  76.             by = bmag * SIN(ellipsearg + 3.14 / 2)
  77.         CASE 19712
  78.             ellipsearg = ellipsearg + .1
  79.             ax = amag * COS(ellipsearg)
  80.             ay = amag * SIN(ellipsearg)
  81.             bx = bmag * COS(ellipsearg + 3.14 / 2)
  82.             by = bmag * SIN(ellipsearg + 3.14 / 2)
  83.     END SELECT
  84.  
  85.     ' Intersections
  86.     a2 = ax ^ 2 + ay ^ 2
  87.     b2 = bx ^ 2 + by ^ 2
  88.     av = ax * vx + ay * vy
  89.     bv = bx * vx + by * vy
  90.     rbx = d - x0
  91.     rby = b - y0
  92.     adbr = ax * rbx + ay * rby
  93.     bdbr = bx * rbx + by * rby
  94.     aa = av ^ 2 / a2 ^ 2 + bv ^ 2 / b2 ^ 2
  95.     bb = 2 * (av * adbr / a2 ^ 2 + bv * bdbr / b2 ^ 2)
  96.     cc = adbr ^ 2 / a2 ^ 2 + bdbr ^ 2 / b2 ^ 2 - 1
  97.     arg = bb ^ 2 - 4 * aa * cc
  98.     IF (arg > 0) THEN
  99.         alpha1 = (-bb + SQR(arg)) / (2 * aa)
  100.         alpha2 = (-bb - SQR(arg)) / (2 * aa)
  101.         x1 = alpha1 * vx + d
  102.         x2 = alpha2 * vx + d
  103.         y1 = alpha1 * vy + b
  104.         y2 = alpha2 * vy + b
  105.     ELSE
  106.         x1 = -999
  107.         y1 = -999
  108.         x2 = -999
  109.         y2 = -999
  110.     END IF
  111.  
  112.     GOSUB draweverything
  113.  
  114.     _LIMIT 60
  115.     _DISPLAY
  116.  
  117.  
  118. draweverything:
  119. PAINT (1, 1), 15
  120. COLOR 0, 15
  121. LOCATE 1, 1: PRINT "LClick=Move ellipse, RClick=Move line, Scroll=Tilt line, Arrows=Shift ellipse"
  122. FOR alpha = -20 TO 20 STEP .001
  123.     x = alpha * vx + d
  124.     y = alpha * vy + b
  125.     CALL ccircle(xorig + x * xzoom, yorig + y * yzoom, 1, 1)
  126. FOR t = 0 TO 6.284 STEP .001
  127.     x = x0 + ax * COS(t) + bx * SIN(t)
  128.     y = y0 + ay * COS(t) + by * SIN(t)
  129.     CALL ccircle(xorig + x * xzoom, yorig + y * yzoom, 1, 4)
  130. CALL ccircle(xorig + x1 * xzoom, yorig + y1 * yzoom, 10, 1)
  131. CALL ccircle(xorig + x2 * xzoom, yorig + y2 * yzoom, 10, 1)
  132.  
  133. SUB cline (x1, y1, x2, y2, col)
  134.     LINE (_WIDTH / 2 + x1, -y1 + _HEIGHT / 2)-(_WIDTH / 2 + x2, -y2 + _HEIGHT / 2), col
  135.  
  136. SUB ccircle (x1, y1, r, col)
  137.     CIRCLE (_WIDTH / 2 + x1, -y1 + _HEIGHT / 2), r, col
« Last Edit: February 01, 2020, 06:47:56 PM by STxAxTIC »
TOXIC

Offline bplus

  • Forum Resident
  • Posts: 6858
  • b = b + ...
Re: Ellipse Intersect Line
« Reply #2 on: February 02, 2020, 12:18:08 AM »
Ah, having done similar, something I can understand and thus appreciate,

A+ man!

Offline STxAxTIC

  • Library Staff
  • Forum Resident
  • Posts: 1062
  • TOXIC
Re: Ellipse Intersect Line
« Reply #3 on: February 02, 2020, 08:06:50 PM »
For completeness, attached is the derivation of the math enclosed above.

To re-derive what's going on in bplus's code next door, let a=b=r and simplify.
« Last Edit: February 02, 2020, 08:30:35 PM by STxAxTIC »
TOXIC