We all love ODF, the best and the most vendor-neutral file-format in the Universe and its surroundings. But for sure, we have some spots where we would prefer it to be somehow less cumbersome. My favourite spot is the need to compute the bounding box of the path element when one writes the `draw:path`

into an OpenDocument Graphics file. Without proper `svg:x`

, `svg:y`

, `svg:height`

, `svg:width`

and `svg:viewBox`

values the path will not be correctly placed.

Computing bounding boxes is not so complicated work when everything is a polygon (which is the case in the internal model of `basegfx`

module), but it becomes a bit more complicated if an application wants to generate paths including elliptical arcs.

I hit this problem about a year ago when I was working during my hackweek on an improvement of `libwpg`

. I tried first to implement the bounding box of an elliptical arc the lazy hacker way, by googling for what other people did. And to my surprise, there is a huge vacuum in what concerns computation of a bounding box of the "`A`

" element of an SVG path. So, I settled for the lazy hacker's plan B: I abandoned the idea by saying I will handle it later, in the first moment, and by implementing a suboptimal solution in the second time. But, since Eilidh did some spectacular progress last week in handling elliptic arcs, my lazyness became the bottle-neck of the progress. So, it was time to remember those nice times when I was warming the benches of the University, dust off my knowledge of mathematical analysis and analytical algebra (or the lack thereof), and try to compute the bounding box of an elliptical arc myself.

And for the purpose of people that might be as lazy as me, I decided to fight my lazyness the second time to give Uncle Google the opportunity to spit out something meaningful, when someone asks it about "Bounding box of an elliptical arc". Here are the notes:

**Compute extremes using parametric description of ellipse**

Let us start from this parametric description of an ellipse:

`x(theta) = cx + rx*cos(theta)*cos(phi) - ry*sin(theta)*sin(phi)`

`y(theta) = cy + rx*cos(theta)*sin(phi) + ry*sin(theta)*cos(phi)`

where `cx`

and `cy`

are the coordinates of the centre of the ellipse, `rx`

and `ry`

are the radii and `phi`

is the x-axis-rotation.

To compute the bounding box of the whole ellipse we need to find for which value of `theta`

the above mentioned functions reach the local extremes. It means where the first derivatives of `x`

and `y`

according to `theta`

are zero. We will get this two equations:

`0 = -rx*sin(theta)*cos(phi) - ry*cos(theta)*sin(phi)`

`0 = -rx*sin(theta)*sin(phi) - ry*cos(theta)*cos(phi)`

which give us two solutions for `x`

:

`theta = -atan(ry*tan(phi)/rx)`

and `theta = M_PI -atan(ry*tan(phi)/rx)`

and two solutions for `y`

:

`theta = atan(ry/(tan(phi)*rx))`

and `theta = M_PI + atan(ry/(tan(phi)*rx))`

**Compute the center of the ellipse**

Since we know now the values of `theta`

describing the extremes of our ellipse, we can compute the `x`

and `y`

values of the bounding box of the whole ellipse. Just to do that, we still need to know the coordinates of the center of the ellipse, `cx`

and `cy`

.

The computation of the center of the ellipse is pretty well described in the Appendix F.6.5 of the SVG standard and does not need to be reproduced here. Just note that for this we need the coordinates of the starting point of the arc that correspond to the end point of the previous path segment.

**Determine the bounding box of the whole ellipse**

Compute the `xmin`

, `xmax`

, `ymin`

and `ymax`

using the values of `theta`

for the local extremes and the newly computed `cx`

and `cy`

coordinates. Like this not only we will know the bounding box of the whole ellipse, but we will also know which value of `theta`

corresponds to maximum and which one to minimum. This knowledge will be later valuable for determining the tightest possible bounding box of a given elliptical arc.

**Tightest possible bounding box**

By calculation of the bounding box of the whole ellipse, we now know the rectangle that will contain the ellipse and thus our elliptical arc too. Nonetheless, this rectangle is too big for our arc. So, the next thing is to trim it down so that it becomes the tightest possible rectangle that will still contain the whole arc.

For this task, we will use the polar coordinates rather then the cartesian ones. The principle is that if any of the points corresponding to `xmin`

, `xmax`

, `ymin`

or `ymax`

of the whole ellipse, lie on the arc they will be be the extremes of the arc too. Nevertheless, if for instance the point `ymin`

does not lie on the arc, the new `ymin`

will be the minimum of the `y`

coordinates of the starting and ending points. In the same way, if the point `xmax`

does not lie on the arc, the new `xmax`

will be the maximum of the `x`

coordinates of the starting and ending points. Whether an extreme does or does not lie on our arc is something trivial to see once the arc is drawn, to determine it programatically will require some efforts.

First, we will compute the coordinates of the points where the whole ellipse touches the bounding box using the parametric description of the ellipse and the values of the `theta`

that we found out in the previous steps. And for determination whether they lie or not on our arc we will use their position in polar coordinates. We will thus need to compute the angles with the x-axis of the lines going through the center of the ellipse and our extreme points. In other terms, we will compute the angle between vector `(1,0)`

and vector `(x`

._{extreme}-cx, y_{extreme}-cy)

The formula for computing the angle between two vectors is known and mentioned *inter alia* as formula F.6.5.4 of the SVG standard. Generally, the expression to calculate the angle between a vector `(ax,ay)`

and a vector `(bx,by)`

is:

`(ax * by > ay * bx ? 1.0 : -1.0) * acos( (ax * bx) + (ay * by) / ( sqrt(ax * ax + ay * ay) * sqrt(bx * bx + by * by) ) )`

.

But since we already know that the first vector is `(1,0)`

, we can simplify it:

`(by > 0.0 ? 1.0 : -1.0) * acos( bx / sqrt(bx * bx + by * by) )`

, which could be eventually simplified to `atan(by / bx)`

, but this expression has a potential division by zero and the code would have to handle those border situations in a special way.

Once we know the angles of the extremes, we still need to calculate the angles of the starting and the end points or our arc using exactly the same formula. So we get `angle1`

corresponding to our starting point and `angle2`

corresponding to the endpoint. It is necessary to normalize all angles so that they lie in the interval of `[0.0, 2.0*M_PI)`

.

In case the `sweep`

flag is 0, the angles are decreasing when the ellipse is drawn. But, for the computation of bounding box the direction of rotation is irrelevant and only complicates the situation. So we swap the angles if the `sweep`

flag is not set. In this way, we can just check for the absence of the extreme points on our elliptical arc, rotating from `angle1`

to `angle2`

. Nevertheless, we have another difficulty with the fact that the angle of 0 radians is the same as the one of 2*M_PI radians. This passage through the 2*M_PI / 0 border is not very easy to handle directly. That is why we swap the points in case where `angle1 > angle2`

and will not look in this case for absence of the extreme points on the arc, but for their presence on the complement arc that would close the ellipse.

And as my teachers used to say: "Grey is the theory, but green is the tree of life," here is what it looks like in a plain C++:

**#include <algorithm>****#include <cmath>****#ifndef M_PI****#define M_PI 3.14159265358979323846****#endif****inline** **double** getAngle(**double** bx, **double** by)**{**

**return** fmod(2*M_PI + (by > 0.0 ? 1.0 : -1.0) * acos( bx / sqrt(bx * bx + by * by) ), 2*M_PI);**}****void** EllpArcBBox(**double** x1, **double** y1,

**double** rx, **double** ry, **double** phi, **bool** largeArc, **bool** sweep, **double** x2, **double** y2,

**double** &xmin, **double** &ymin, **double** &xmax, **double** &ymax)**{**

**if** (rx < 0.0)

rx *= -1.0;

**if** (ry < 0.0)

ry *= -1.0;

**if** (rx == 0.0 || ry == 0.0) **{**

xmin = (x1 < x2 ? x1 : x2);

xmax = (x1 > x2 ? x1 : x2);

ymin = (y1 < y2 ? y1 : y2);

ymax = (y1 > y2 ? y1 : y2);

**return**;

**}**

**const** **double** x1prime = cos(phi)*(x1 - x2)/2 + sin(phi)*(y1 - y2)/2;

**const** **double** y1prime = -sin(phi)*(x1 - x2)/2 + cos(phi)*(y1 - y2)/2;

**double** radicant = (rx*rx*ry*ry - rx*rx*y1prime*y1prime - ry*ry*x1prime*x1prime);

radicant /= (rx*rx*y1prime*y1prime + ry*ry*x1prime*x1prime);

**double** cxprime = 0.0;

**double** cyprime = 0.0;

**if** (radicant < 0.0) **{**

**double** ratio = rx/ry;

**double** radicant = y1prime*y1prime + x1prime*x1prime/(ratio*ratio);

**if** (radicant < 0.0) **{**

xmin = (x1 < x2 ? x1 : x2);

xmax = (x1 > x2 ? x1 : x2);

ymin = (y1 < y2 ? y1 : y2);

ymax = (y1 > y2 ? y1 : y2);

**return**;

**}**

ry=sqrt(radicant);

rx=ratio*ry;

**}** **else** **{**

**double** factor = (largeArc==sweep ? -1.0 : 1.0)*sqrt(radicant);

cxprime = factor*rx*y1prime/ry;

cyprime = -factor*ry*x1prime/rx;

**}**

**double** cx = cxprime*cos(phi) - cyprime*sin(phi) + (x1 + x2)/2;

**double** cy = cxprime*sin(phi) + cyprime*cos(phi) + (y1 + y2)/2;

**double** txmin, txmax, tymin, tymax;

**if** (phi == 0 || phi == M_PI) **{**

xmin = cx - rx;

txmin = getAngle(-rx, 0);

xmax = cx + rx;

txmax = getAngle(rx, 0);

ymin = cy - ry;

tymin = getAngle(0, -ry);

ymax = cy + ry;

tymax = getAngle(0, ry);

**}** **else** **if** (phi == M_PI / 2.0 || phi == 3.0*M_PI/2.0) **{**

xmin = cx - ry;

txmin = getAngle(-ry, 0);

xmax = cx + ry;

txmax = getAngle(ry, 0);

ymin = cy - rx;

tymin = getAngle(0, -rx);

ymax = cy + rx;

tymax = getAngle(0, rx);

**}** **else** **{**

txmin = -atan(ry*tan(phi)/rx);

txmax = M_PI - atan (ry*tan(phi)/rx);

xmin = cx + rx*cos(txmin)*cos(phi) - ry*sin(txmin)*sin(phi);

xmax = cx + rx*cos(txmax)*cos(phi) - ry*sin(txmax)*sin(phi);

**if** (xmin > xmax) **{**

std::swap(xmin,xmax);

std::swap(txmin,txmax);

**}**

**double** tmpY = cy + rx*cos(txmin)*sin(phi) + ry*sin(txmin)*cos(phi);

txmin = getAngle(xmin - cx, tmpY - cy);

tmpY = cy + rx*cos(txmax)*sin(phi) + ry*sin(txmax)*cos(phi);

txmax = getAngle(xmax - cx, tmpY - cy);

tymin = atan(ry/(tan(phi)*rx));

tymax = atan(ry/(tan(phi)*rx))+M_PI;

ymin = cy + rx*cos(tymin)*sin(phi) + ry*sin(tymin)*cos(phi);

ymax = cy + rx*cos(tymax)*sin(phi) + ry*sin(tymax)*cos(phi);

**if** (ymin > ymax) **{**

std::swap(ymin,ymax);

std::swap(tymin,tymax);

**}**

**double** tmpX = cx + rx*cos(tymin)*cos(phi) - ry*sin(tymin)*sin(phi);

tymin = getAngle(tmpX - cx, ymin - cy);

tmpX = cx + rx*cos(tymax)*cos(phi) - ry*sin(tymax)*sin(phi);

tymax = getAngle(tmpX - cx, ymax - cy);

**}**

**double** angle1 = getAngle(x1 - cx, y1 - cy);

**double** angle2 = getAngle(x2 - cx, y2 - cy);

**if** (!sweep)

std::swap(angle1, angle2);

**bool** otherArc = **false**;

**if** (angle1 > angle2) **{**

std::swap(angle1, angle2);

otherArc = **true**;

**}**

**if** ((!otherArc && (angle1 > txmin || angle2 < txmin)) || (otherArc && !(angle1 > txmin || angle2 < txmin)))

xmin = x1 < x2 ? x1 : x2;

**if** ((!otherArc && (angle1 > txmax || angle2 < txmax)) || (otherArc && !(angle1 > txmax || angle2 < txmax)))

xmax = x1 > x2 ? x1 : x2;

**if** ((!otherArc && (angle1 > tymin || angle2 < tymin)) || (otherArc && !(angle1 > tymin || angle2 < tymin)))

ymin = y1 < y2 ? y1 : y2;

**if** ((!otherArc && (angle1 > tymax || angle2 < tymax)) || (otherArc && !(angle1 > tymax || angle2 < tymax)))

ymax = y1 > y2 ? y1 : y2;**}**