# Determining which octant has a specific point

I am interested in octants in plane geometry, which basically involves dividing a plane into 8 equal parts according to a certain point.

Octants are particularly useful while implementing the Bresenham's algorithm.

## Problem

Given two points A = (x1, y1) and B = (x2, y2), how to determine in which octant belongs B according to A?

Counting octants is performed as shown: starting from the upper right and iterating anticlockwise.

## Naive implementation

def get_octant(pA, pB):
AX = pA[0]
AY = pA[1]
BX = pB[0]
BY = pB[1]

if BX > AX:
if BY > AY:
if (BY - AY) < (BX - AX):
octant = 0
else:
octant = 1
else:
if (AY - BY) < (BX - AX):
octant = 7
else:
octant = 6
else:
if BY > AY:
if (BY - AY) < (AX - BX):
octant = 3
else:
octant = 2
else:
if (AY - BY) < (AX - BX):
octant = 4
else:
octant = 5

return octant


Would it not possible to find a cleverer implementation?

Using the symmetry properties, I am convinced that it is possible to make a shorter and more readable code. Unfortunately, I do not find it myself, so does anyone have an idea, please?

• Unrelated: to support A.x syntax, you could use collections.namedtuple – jfs Jul 1 '15 at 12:28
• A small suggestion: rather than assigning Ax, Ay, ... on individual lines, consider using tuples unpacking: Ax, Ay = pA; Bx, By = pB. – alexwlchan Jul 2 '15 at 12:07

I am convinced that it is possible to make a shorter and more readable code.

To implement the discrete analog of atan2((B.y - A.y), (B.x - A.x)):

x, y  = (B.x - A.x), (B.y - A.y)
octant = [([1, 2], [8, 7]), ([4, 3], [5, 6])][x < 0][y < 0][abs(x) < abs(y)]


It is shorter, I don't know about readability. It is easy to write tests that check all possible variants including the behavior on the boundaries.

Yes, there are exploitable symmetries.

def get_octant(pA, pB):
dx, dy = pB[0] - pA[0], pB[1] - pA[1]
octant = 0
if dy < 0:
dx, dy = -dx, -dy  # rotate by 180 degrees
octant += 4
if dx < 0:
dx, dy = dy, -dx  # rotate clockwise by 90 degrees
octant += 2
if dx < dy:
# no need to rotate now
octant += 1
return octant

• (This doesn't have the same behavior as your code on octant boundaries. It's not clear whether this is an issue.) – David Eisenstat Jun 30 '15 at 17:58
• you beat me to it, as often happens. =) However my solution was slightly less obvious as to why it worked, and slightly more code. +1 – user2566092 Jun 30 '15 at 18:04

I would go with finding the angle between AB and a line parallel to the oX axis. You can simulate the line parallel to oX that passes through A with AC where C(A.x+k,A.y) (for any k>0). Now, since ABC forms a triangle, you can exploit the cosine law, hence the angle <BAC being acos((AB^2 + AC^2 - BC^2)/(2*AB*AC)) this will give you an angle in the range [0,pi]. In order to adjust it, you would like to know whether B is above or below AC. This can be computed with sign(B.y - A.y) which gives 1 if B is above and -1 otherwise.

Putting it all together: ((-180*(sign(B.y - A.y)-1) + acos((AB^2 + AC^2 - BC^2)/(2*AB*AC)) * sign(B.y - A.y)) div 45)+1. This will give exactly the values from your image.

However, this doesn't mean that this one liner will execute faster then the nested if/else you have, reason for this: implementation of acos() is non trivial; no optimisations e.g. branch prediction.

• This is the math I was looking for, thank you! I accepted Sebastian's answer as it is more understandable and probably faster, but yours is what I had in mind while publishing this question. – Delgan Jul 3 '15 at 16:09

You can look at the computation of the octant number (q) of a paper that I co-authored:

import numpy as np