Draw largest rectangle within a limited shape

[hs800150]

Ok so this lisp might be too complicated but I would like to know if anyone knows a lisp to automate the largest rectangle within a space. I have attached a picture of what I mean. The red rectangle has to be a rectangle (90 degrees between each angle) and stay within the boundaries of the white shape. Lets say I can choose a width say 40' and the lisp automatically finds the largest length based on trial and error. Or if there is a lisp to find both the largest width and length that would work too. 

[john.uhden]

In your "Lot Fit1,png" it appears that the concave curve represents the front setback on a cul-de-sac.  It also appears that the largest rectangle might be achieved by setting the front building line back even farther because the lot width increases toward the rear, which of course makes the programming more complex.

[hak_vz]

@hs800150  Have you tried provided solutions?

[jtoverka]

I believe the solution to this problem is quite simple based upon my understanding of the constraints.

Constraints:

1. Must be the maximum area rectangle.

2. The rectangle must be within the white polyline.

These constraints are quite simple.

Now comes the assumptions on the problem statement. Hopefully my assumptions are valid and correct. After looking at the examples provided, it appears that at a minimum, two opposing corners are always touching the polyline. Additionally, it appears that the maximum rectangle always covers the centroid of the polyline.

Assumptions

1. Two opposing corners are always touching the polyline.

2. No padding has been applied to the polyline (this wouldn't complicate much, ignored for simplicity).

3. The rectangle encompasses the centroid.

Code for this can be found here.

Lee Mac's polygon centroid 

Lee Mac's geometric functions 

A central limitation to the polygon centroid is that it requires straight lines for the polyline. I believe this can be modified to accommodate other types of segments like arcs.

If you have a polyline with arc segments the following code can be used to create information on this polyline:

Lee Mac's polyinfo 

It has the code you will need to create a list of points out of the polyline.

Now, for the algorithm:

1. Convert the polyline into a series of points in a list.

   a) The points are generated based on a precision parameter e.g. 5000 points

2. Use another precision parameter for the rotation of the rectangle e.g. 360 (1 degree precision)

3. Three nested loops

   a) This is not required as it will create duplicate calculations, and you can create a workaround for a reduction in the big O. The nested point loops ignoring the rotation loop is O(n^2). Where as a clever form of these loops is O(n).

(foreach point1 pointList
  (foreach point2 pointList
    (repeat 360
      ; process here
    )
  )
)

4. Calculate the rectangle area, check the following:

    a) is the calculated rectangle within the white line? use the centroid as a reference point. This may require a function that does more loops to verify if a line on the rectangle does not intersect the polyline, but does not count the end points of the rectangle onto the polyline.

    b) If the rectangle area is larger than the maximum found area, then store the points in variables.

After this process, you're done!

You should have the maximum found area within the precision parameters you have listed.

[doaiena]

You should not make assumptions, based on just a few examples. You should try to make your code as generic as possible, to cover all possible situations you can think of. I have 2 examples, that make your assumptions invalid.

As for the solution. Me and @hak_vz have posted only partial solutions at best, with loads of limitations. I don't think it's as simple as you have stated. Maybe you should post some code and prove me wrong?

[jtoverka]

Assumptions are important, however, you are correct. These assumptions are invalid. I will hold onto one of these assumptions, that is a closed polyline.

The best we can do at this point is to use ray tracing. A single scan can yield a list of intersections. Do another scan at 90 degrees at these intersection points. To figure out what range of intersection points are inside or outside the polyline, refer to the figure below. Additionally, further scans at differing angles are required to find the maximum area rectangle.

The odd number intersection point indicated inside the polygon, whereas the even number intersection point is outside the polygon. This is part of the point in polygon problem. 

You can use a literal ray object and the intersectWith method on the polyline. Move the basepoint of the ray to scan, and use the direction vector for the angle of the scan. The distance that you move the ray will be your precision.

Although using a literal moving ray object is less efficient than doing the actual calculation, it is the simplest solution.

[Kent1Cooper]

---

@jtoverka wrote:

....

The odd number intersection point indicated inside the polygon, whereas the even number intersection point is outside the polygon. ....

---

That is another not-always-valid assumption.  It requires that the ray does not exactly touch  any vertex on its way through:

Here, the point at left is outside, but the ray finds an odd  number of intersections.

[marko_ribar]

---

@Kent1Cooper wrote:

---

@jtoverka wrote:

....

The odd number intersection point indicated inside the polygon, whereas the even number intersection point is outside the polygon. ....

---

That is another not-always-valid assumption.  It requires that the ray does not exactly touch  any vertex on its way through:

Here, the point at left is outside, but the ray finds an odd  number of intersections.

---

Gilles sub function may be used to overcome this issue you pointed Kent... Look closely end portion of this code :

; Lee Mac Point Inside Curve
(defun LM:Inside-p ( pt ent / unit v^v _GroupByNum fd1 fd2 par lst nrm obj tmp )

  (vl-load-com)

  (defun *error* ( errmsg )
    (if (not (wcmatch errmsg "Function cancelled,quit / exit abort,console break"))
      (princ (strcat "\nError: " errmsg))
    )
    (vla-put-color obj acYellow)
    (princ)
  )

  (defun unit ( v / d )
    (if (not (equal (setq d (distance '(0.0 0.0 0.0) v)) 0.0 1e-8))
      (mapcar '(lambda ( x ) (/ x d)) v)
    )
  )

  (defun v^v ( u v )
    (list
      (- (* (cadr u) (caddr v)) (* (caddr u) (cadr v)))
      (- (* (caddr u) (car v)) (* (car u) (caddr v)))
      (- (* (car u) (cadr v)) (* (cadr u) (car v)))
    )
  )

  (defun _GroupByNum ( l n / r )
    (if l
      (cons (reverse (repeat n (setq r (cons (car l) r) l (cdr l)) r)) (_GroupByNum l n))
    )
  )

  (if (= (type ent) 'VLA-OBJECT)
    (setq obj ent
          ent (vlax-vla-object->ename ent))
    (setq obj (vlax-ename->vla-object ent))
  )

  (if (vlax-curve-isplanar ent)
    (progn
      (setq fd1 (vlax-curve-getfirstderiv ent (setq par (vlax-curve-getstartparam ent))))
      (while (equal fd1 (setq fd2 (vlax-curve-getfirstderiv ent (setq par (+ par 0.001)))) 1e-3))
      (setq nrm (unit (v^v fd1 fd2)))
      (setq lst
        (_GroupByNum
          (vlax-invoke
            (setq tmp
              (vlax-ename->vla-object
                (entmakex
                  (list
                    (cons 0 "RAY")
                    (cons 100 "AcDbEntity")
                    (cons 100 "AcDbRay")
                    (cons 10 pt)
                    (cons 11 (trans '(1. 0. 0.) nrm 0))
                  )
                )
              )
            )
            'IntersectWith obj acextendnone
          ) 3
        )
      )
      (vla-delete tmp)
      ;; gile:
      (and
        lst
        (not (vlax-curve-getparamatpoint ent pt))
        (= 1 (rem (length (vl-remove-if (function (lambda ( p / pa p- p+ p0 )
                                                    (setq pa (vlax-curve-getparamatpoint ent p))
                                                    (and (setq p- (cond ((setq p- (vlax-curve-getPointatParam ent (- pa 1e-8)))
                                                                         (trans p- 0 nrm)
                                                                        )
                                                                        ((trans (vlax-curve-getPointatParam ent (- (vlax-curve-getEndParam ent) 1e-8)) 0 nrm)
                                                                        )
                                                                  )
                                                         )
                                                         (setq p+ (cond ((setq p+ (vlax-curve-getPointatParam ent (+ pa 1e-8)))
                                                                         (trans p+ 0 nrm)
                                                                        )
                                                                        ((trans (vlax-curve-getPointatParam ent (+ (vlax-curve-getStartParam ent) 1e-8)) 0 nrm)
                                                                        )
                                                                  )
                                                         )
                                                         (setq p0 (trans pt 0 nrm))
                                                         (<= 0 (* (sin (angle p0 p-)) (sin (angle p0 p+)))) ;; LM Mod
                                                    )
                                                  )
                                        ) lst
                          )
                  ) 2
             )
        )
      )
    )
    (prompt "\nReference curve isn't planar...")
  )
)

[jtoverka]

You can either increase the precision, or add optimization to account for the room between the precision points. The simplest solution would be to increase the precision. I really don't think there is need to fret over a couple of square inches out of what? acres?

[john.uhden]

We've been through this before. You can't assume that the odd/even
approach will give you the results you expect. A ray cast from a point
outside the polygon might strike exactly one vertex. And a ray cast from a
point inside a polygon might intersect once but also run through exactly
one vertex.
I recently posted a solution. I think I called it "@inside revisited."

[jtoverka]

Why would a ray cast be inside the polyline? You would take the bounding
box prior casting Ray's. You would obviously know that a ray is cast
outside a polyline. And so what if it lands exactly on a vertex? The next
ray that is like a micron away will be clear from the prior vertex.

Not only that, attacking it with differing angles of ray's makes it very
unlikely that a solution could not be found.

[hak_vz]

As it seams to me OP has either gone with solutions posted before ( @CodeDing @doaiena and @hak_vz) or he has simply forgotten about this post.

It is no matter what method you want to apply but to come out with a code that is simple and that works. Try to write down some code.

For the purpose OP has asked featured codes are more than satisfiying. Each of them needs final touches.

I would in my code include entmake instead command, include testing for sharp vertexes since iterative code can skip it (explained by @Kent1Cooper in answer to your post) and that is all.

[john.uhden]

I am not quoting you exactly, but I think you were saying that you could
determine whether a point was inside or outside of a polygon by ray
casting, with one intersection meaning inside and none or two meaning
outside. I was just pointing out that that presumption is not 100%
correct. Plus I failed to mention that an intersection from a point
outside could just be a tangent point on a bulged segment.
But if you have all that covered in the rest of your code, then please
accept my apology for interrupting your progress.

[CodeDing]

I've been working on this in my free time and it's a great project.

Attached is a short clip of my progress.

I'm taking a sort-of pixelated approach. Depending on function inputs the precision can grow, but in essence, a grid of squares is created around the object at different angles, then we can analyze which squares are inside the polygon and ultimately determine which 'area' of combined squares provides the largest area.

I had it draw the squares and pause each time so that we can see the current progress.

Next steps for me are:

1) Determining which squares are inside the polygon (I'll ultimately use the Sum of Angles method)

2) Using only squares inside of polygon to determine largest area.

My goal with this method is to get CLOSE to the largest shape, not exact, because as we can all see it's a very hard thing to determine.

Video attached. (I don't have ScreenCast Lol)

Best,

~DD

[john.uhden]

The Sum of Angles approach is the fastest. It should work just fine for
polygon lots. Just remember to use a parameter increment of 0.5 to account
for bulges.

[CodeDing]

Ok, I'm down to my last item on my checklist. It's just difficult to grasp, so I'm asking if anybody else has a function like this or can help. I've been able to successfully (mostly) determine which 'pixels' are inside of my polygon.

Example:

 

So, last step is this...

Given a 2D boolean matrix (in our case a list of lists (with equal length)), return the largest rectangle of T values.

Here is a video describing the implementation (with binary):

https://www.youtube.com/watch?v=g8bSdXCG-lA

And here is some code in other languages:

https://www.geeksforgeeks.org/maximum-size-rectangle-binary-sub-matrix-1s/

However, since this 'area' needs to be correlated to another list (or matrix), I also need the 'coordinates' returned with this function. And I cannot grasp how begin writing it.

Here are two examples of expected input / output:

---IN---
(setq matrix '(
  (nil  T   T  nil nil)
  (nil  T   T  nil nil)
  (nil  T   T   T  nil)
  (nil nil  T  nil nil)
  (nil nil nil  T  nil)
))
(LargestArea matrix 0) ;<-- '0' is the largest previously found area,
                       ;--- since this function is inside a loop.
---OUT---
(6 (2 1) (0 2))
;;6 is largest area
;;(2 1) is Lower Left 'coordinate'
;;(0 2) is Upper Right 'coordinate'

---IN---
(setq matrix '(
  (nil nil  T nil nil)
  (nil nil  T  nil nil)
  (nil  T   T   T   T)
  ( T   T   T   T   T)
  (nil nil nil  T  nil)
))
(LargestArea matrix 6) ;<-- '6' is the largest previously found area,
                       ;--- since this function is inside a loop.
---OUT---
(8 (3 1) (2 4))
;;8 is largest area
;;(3 1) is Lower Left 'coordinate'
;;(2 4) is Upper Right 'coordinate'

Any help would be appreciated, but if not I'll still keep chipping away.

Best,

~DD

[john.uhden]

That LL / UR method may be better than what my elder mind was thinking,
being the four corners as 3D or even 2D points. But it might get confusing
if you iterate through various rotations since the LL will ultimately
become the UR, but comparing the two values should take care of that.

[hs800150]

@hak_vz 

I am liking this code so far and works pretty well. I was wondering if this code could be modified some. I want to be able to select multiple polygons and then I could let the computer calculate each one (one by one) based on the selections I make. For example:

1. The prompt asks: "Select Polygon" Action: I select the polygon -> Then

2. The prompt asks: "Select a point on a polygon side that rectangle edge is aligned to" Action: I select the edge -> Then

2. The prompt asks: "Would you like to select another polygon? [y/n]" Action: y -> Then

3. 1. and 2. repeat for another polygon

4. The prompt asks: "Would you like to select another polygon? [y/n]" Action: n -> Then computer performs the box calculations for 2 polygons.

Do you think this is possible?

Thanks

[hak_vz]

I'll try to modify a script to perform that way. I'm working on another project and Easter holidays are here.

I don't see a reason to do that since time of execution will be similar.

If you want me to create batch process that would perform action on let say 50 parcels we can do as follow. For each side on a polygon that has to be tested you can place a circle at midpoint (in a separate layer so you later can erase them without to much fuss). For each circles code will search  for a polygon that it intersects with and create max rectangle.  

Circle can be also placed inside a polygon near midpoint of a side we want to place our max rectangle, to avoid intersecting with other polygons.

[john.uhden]

Of course it's possible...

(while (not done)...
)

[hs800150]

@hak_vz 

How about instead of a circle I draw one polyline through all of the boxes and whatever the polyline intersects with will draw the largest box?