If you are just starting with Python I hope the small problems in this blog are helpful. All you need to solve or follow them is very basic knowledge of Python and curiosity.
Each post uses trinkets to run python on the browser. You can modify and run the code on the blog post itself.
Showing posts with label geometry. Show all posts
Showing posts with label geometry. Show all posts

Saturday, April 18, 2015

How many carrots on a line?

You live a simple and happy life in your farm. 
You grow all sorts of vegetables in the fields.
Some days of the year, you toil the soil all day planting carrots, tomatoes, lettuce, cabbage, beet, pepper and many others. 
At night you are so tired that you sleep deeply and without any dreams.
The rains come and the fields become a deep luxuriant green with spots of red, yellow, orange, purple and patches of mint green.
One day while harvesting that season's carrots, you sit midday on a rock for a brief respite. Lines and lines of carrots extend in every direction. 
The sun shines strongly above. You are thirsty and tired. Suddenly you see some leaves move slightly. As in a dream (perhaps you did indeed take a nap?) you see a neatly dressed carrot man approach you.

"Good day sir" says the carrot man. After a short pause as if expecting you to say something he continues "we see you working all day pulling our neighbors from the soil and we ask each other how long before we too will be taken from our homes?".

"I have never thought of it that way. I don't even know how many of you there are?" you mumble a bit embarrassed.

"There is a very simple way to estimate that sir. Imagine that there is one carrot being living at each integer set of coordinates (x,y) of this field.. Then you certainly can figure out how many carrots lie in the line between any two carrots at coordinates (x1,y1) and (x2,y2)?  If you then extend your reasoning to find how many carrots are in a rectangle that has its diagonal given by the same line segment you can find out how many of us there are."
Having completed this discourse the carrot man turned and slowly walked to the neighboring row of carrots, past it and soon he had disappeared among the greenery.
After this conversation you hesitate going back to your carrot pulling task. You try to solve the problem posed by the carrot man instead. How many carrots there are in the line segment between two carrots? And how many carrots lie in a rectangle whose diagonal is given by the same line segment?

How many carrots on a line?

Cross Product

Given two carrots at (x1,y1) and (x2,y2) a simple strategy is the following:
Consider all carrots that lie inside (but not on the border) of a rectangle where the lower left corner is at coordinates min(x1,x2), min(y1,y2) and the top right corner is at coordinates max(x1,x2),max(y1,y2). 
We can check for each of these carrots whether they lie or not in the segment connecting (x1,y1) to (x2,y2). How? 
Form the vector that joins the two ends of the segment. Its components will be (max(x1,x2)-min(x1,x2), max(y1,y2)-min(y1,y2)). 
Then form the vector connecting the lower left corner of the rectangle at (min(x1,x2),min(y1,y2)) to the point inside the rectangle at integer coordinates (i,j). 
This vector has coordinates (i-min(x1,x2)),j-min(y1,y2)). 
The point lies on the segment when the angle between these 2 vectors is zero.

The cross product of these two vectors is the product of their magnitudes times the sin of the angle between them. Therefore the cross product will vanish if these two vectors are collinear, i.e. if the point lies on the segment.
The cross product is a vector parallel to the z axis (the axis perpendicular to the field where the carrots lay) with value

(i-min(x1,x2))*(max(y1,y2)-min(y1,y2)) - (j-min(y1,y2))*(max(x1,x2)-min(x1,x2))


All points inside the rectangle for which this difference is zero, are on the line segment.

Greatest Common Divisor (GCD)

Since we are dealing with integers,  there is another way we can find the number of points using number theory.

Without loss of generality we can always put the lower left corner of the rectangle at the origin of the coordinates (0,0) and the top right corner at the point (x2-x1,y2-y1) (assuming x2>x1 and y2>y1). The rectangle has sides of length (x2-x1) and (y2-y1).


All integer points on this segment must obey the line equation y=(y2-y1)/(x2-x1)x. Or y/x=(y2-y1)/(x2-x1).


Lets suppose that the GCD of (y2-y1) and (x2-x1) is ...GCD. We can write


y = GCD*p/(GCD*q)x

where p and q are relatively prime (coprime), i.e. they have no common divisors.


It is clear then that the following integer points are on the segment: (q,p), (2q,2p), (3q,3p),..., ((GCD-1)x, (GCD-1)y). Therefore there are GCD-1 integer points/carrots in the segment (excluding the endpoints).


For an example consider (x2-x1)=10 and (y2-y1)=15. Then y/x=3*5/(2*5). Pairs (2,3), (4,6), (6,9), (8,12) are all in the segment. The GCD is 5 and there are 4 points on the segment.


This also shows that if  (x2-x1) and (y2-y1) are coprimes (i.e. they dont have any common divisors), there can not be any integer point lying on the segment (other than the endpoints).

Put another way dividing the rectangle of size (x2-x1,y2-y1) into the greatest number of identical rectangles is equivalent to finding the greatest common divisor of (x2-x1) and (y2-y1).

How many carrots in a rectangle?

Lets now look at the second question:
How many carrots are inside the rectangle whose diagonal is the segment from (x1,y1) to (x2,y2)?
The easiest way to count those is to count them as we loop through every point of integer coordinates on or inside the rectangle.

But since the problem is so easy we can try to be clever and come up with a formula for this number.

To find out how many carrots are right on the rectangle border we could use the cross product or gcd methods above applied to each of the four sides of the rectangle.  Actually the top and bottom sides will have the same number of carrots as will the left and right sides so we can just use two sides and multiply by two.
But it is easier to do things manually:
The bottom side will have (x2-x1)+1 points (because of the origin). Adding the top side this gives
2[x2-x1+1]. The left and right side contribute 2[y2-y1+1]. But 4 of the points are counted twice, so
the total number of carrots on the boundaries of the rectangle is
2[x2-x1+1+y2-y1+1-2]=2(x2-x1+y2-y1)

How many carrots will be inside the rectangle? Divide the rectangle by vertical lines parallel to its left and right sides and starting at integer values of x. There will be (x2-x1)-1 of those lines (excluding the left and right sides). For each of those lines there will (y2-y1)-1 points with integer y coordinates. So there will be a total of [x2-x1-1][y2-y1-1] interior carrots.


The total number of carrots in the rectangle (interior and boundary) is then


(x2-x1)(y2-y1)+ (x2-x1+y2-y1)+1 


The following code implements both the cross product and gcd solutions. In particular notice the implementation of gcd using the Euclidean algorithm. The code also counts the number of points inside and on the rectangle and compares it with the formula above.


Monday, March 23, 2015

Intersection of two line segments in the plane

We are given two line segments (in the plane) through their two end point (x,y) coordinates. To simplify the problem lets assume the two segments are either horizontal or vertical and their endpoint coordinates are integers.

The two segments might not intersect at all, they might intersect at a single point or they might have a whole segment in common. The following figure shows those three possibilities. The last pair of segments (red and cyan) overlap in a segment, the pair before the red-cyan pair intersects at a single point and the first three pairs of segments do not intersect. Notice that the overlap segment will either be horizontal or vertical since the two segments are both horizontal or both vertical when a segment overlap occurs.

How would we write a simple python program that returns a string 'NO INTERSECTION', 'POINT INTERSECTION' or 'SEGMENT INTERSECTION' for each of the three possible cases outlined above and returns the intersection segment?

 Lets suppose the endpoint coordinates of the first segment are given as a python list s1:[x1 y1 x2 y2] and similarly for the second segment s2:[x1 y1 x2 y2].

 It is tempting to use a few if statements to check for the different configurations of the two line segments but boundary cases make such an approach hard to get right: Imagine one of the line segments is just a single point in the other segment. In this case we want to return "POINT INTERSECTION" not "SEGMENT INTERSECTION". Yet another case is if the two line segments are both points.

 An easier and better solution is to write code that naturally handles all possible cases. We can first easily determine the overlapping region of the two line segments, and then see if it is a point, a segment, or if there is no overlap.

The x coordinate of the leftmost endpoint of the intersection segment will be

  left = max(min(s1.x1,s1.x2),min(s2.x1,s2.x2)).

That is we search among the endpoints of s1 and s2 for the the one with the largest left coordinate  (the min gets the endpoint in each segment that is the leftmost).

Similarly the x coordinate of the rightmost endpoint of the intersection segment will be

right = min(max(s1.x1,s1.x2),max(s2.x1,s2.x2)).

That is we search among the endpoints of s1 and s2 for the one with the smallest right coordinate (the max gets the endpoint in each segment that is the rightmost).

We can extend this logic to find the y coordinates of the endpoints of the intersection segment:

bottom = max(min(s1.y1,s1.y2),min(s2.y1,s2.y2))

top = min(max(s1.y1,s1.y2),max(s2.y1,s2.y2))

With the coordinates of the segment endpoints we can promptly find out if there is intersection or not. For example if left > right or bottom>top there is no intersection whatsoever. If left=right and bottom=top the 2 segments meet at a single point. Otherwise they overlap in a segment.


Here it is this solution (you can change the code and try different approaches: it will be executed roughly every minute, and you can change this interval by editing the value 60000 to the value most convenient to you)
And a more concise way of coding it up is