# geometry.frink

``` // Formulae for simple plane geometry. // Find the length of the hypotenuse of a right triangle given sides a and b. hypotenuseLength[a, b] := sqrt[a^2 + b^2] // Calculate area of arbitrary triangle given 3 sides of known length. // Calculate using Heron's formula. triangleArea[a, b, c] := {    // Semiperimeter s    s = 1/2 (a + b + c)    area = sqrt[s (s-a) (s-b) (s-c)]    return area } // Find the distance between a point (x0,y0) and a line specified by // two points (x1, y1), (x2, y2) pointToLineDistance[x0, y0, x1, y1, x2, y2] := {    abs[ (x2-x1)(y1-y0) - (x1-x0)(y2-y1)] / sqrt[(x2-x1)^2+(y2-y1)^2 ] } // Returns true if a point is inside a polygon. // Adapted from // http://www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/pnpoly.html // // points is an array of two-dimensional arrays of points. // x,y is the point to test. // // Returns: //   true if [x,y] is inside the polygon, false otherwise. pointInPolygon[points, x, y] := {    c = false    i = 0    npol = length[points]    j = npol-1    while i < npol    {       xi = points@i@0       xj = points@j@0       yi = points@i@1       yj = points@j@1       if ((((yi<=y) and (y<yj)) or ((yj<=y) and (y<yi))) and (x < (xj - xi) * (y - yi) / (yj - yi) + xi))          c = !c       j = i       i = i + 1    }    return c; } // Calculates the area of a two-dimensional polygon. // // points is an array of [x,y] coordinates // // Returns: [cx, cy] // //  See:  http://en.wikipedia.org/wiki/Polygon#Area_and_centroid // //  The Frink polygon classes now implement area and centroid calculations // using this code. polygonArea[points] := {    size = length[points]    selfClosed = (points@0 == points@(size-1))    if selfClosed       n = size - 1    else       n = size    if size == 0       return 0    if size <= 2       return 0 * points@0@0 * points@0@1    // Maintain units of measure.    xs = 0 * points@0@0    ys = 0 * points@0@1    ps = xs ys    xs = xs * ps    ys = ys * ps        for i = 0 to n-1    {       xi = points@i@0       yi = points@i@1       ii = i+1       if ii == n and !selfClosed  // Make last point first point          ii = 0              xi1 = points@(ii)@0       yi1 = points@(ii)@1       pp = (xi yi1 - xi1 yi)       ps = ps + pp    }    A = 1/2 ps    return abs[A] } // Calculates the centroid of a two-dimensional polygon. // // points is an array of [x,y] coordinates // // Returns: [cx, cy] // //  See:  http://en.wikipedia.org/wiki/Polygon#Area_and_centroid polygonCentroid[points] := {    size = length[points]    selfClosed = (points@0 == points@(size-1))    if selfClosed       n = size - 1    else       n = size    if size == 0       return undef    if size == 1       return points@0    // Maintain units of measure.    xs = 0 * points@0@0    ys = 0 * points@0@1    ps = xs ys    xs = xs * ps    ys = ys * ps        for i = 0 to n-1    {       xi = points@i@0       yi = points@i@1       ii = i+1       if ii == n and !selfClosed  // Make last point first point          ii = 0              xi1 = points@(ii)@0       yi1 = points@(ii)@1       pp = (xi yi1 - xi1 yi)       xs = xs + (xi + xi1) pp       ys = ys + (yi + yi1) pp       ps = ps + pp    }    A = 1/2 ps    cx = xs/(6 A)    cy = ys/(6 A)    return[cx, cy] } // Rotate the point [x,y] around the point [cx,cy] by the specified angle // and return the coordinates of the new point, [x2, y2] rotateAroundPoint[x,y, cx, cy, angle] := {    ca = cos[angle]    sa = sin[-angle]    dx = x-cx    dy = y-cy    return [cx + ca dx - sa dy, cy + sa dx + ca dy] } // Finds the intersection[s] between a line and a circle given two points on // the line [x1, y1], [x2, y2], the center of the circle [cx, cy], // and the radius of the circle. // //  Returns: //   an array, each element containing an array with [x,y] coordinates for //  each intersection point.  There can  //  be 0, 1, or 2 points of intersection.  For example, this may return //  something like  [[1,1], [2,1]] // // See: // http://mathworld.wolfram.com/Circle-LineIntersection.html circleLineIntersections[x1, y1, x2, y2, cx, cy, r] := {    // First shift coordinates so the center of the circle is at 0,0    x1a = x1 - cx    y1a = y1 - cy    x2a = x2 - cx    y2a = y2 - cy    dx = x2 - x1    dy = y2 - y1    dr2 = dx^2 + dy^2    D = x1a y2a - x2a y1a    delta = r^2 dr2 - D^2    // If delta < 0 there is no intersection.    if delta < 0       return new array    sgndy = dy < 0 ? -1 : 1    sqrtDelta = sqrt[delta]    Ddy = D dy    Ddx = - D dx        // if delta == 0 there is only one point of intersecion    if delta == 0       return [[Ddy / dr2 + cx, Ddx / dr2 + cy]]        rightx = sgndy dx sqrtDelta    righty = abs[dy] sqrtDelta    return [[(Ddy + rightx) / dr2 + cx, (Ddx + righty) / dr2 + cy],            [(Ddy - rightx) / dr2 + cx, (Ddx - righty) / dr2 + cy]] } /** This performs an arbitrary affine transformation on the point [x,y]     using the affine coordinates a,b,c,d,e,f.     where the transformation matrix is defined by:     x'    a c e     x     y' =  b d f  *  y     1     0 0 1     1     (This is the order of the coordinates in the SVG and HTML5 standards.)     which gives the equation:     x' = a x + c y + e     y' = b x + d y + f     returns:  [x', y'] */ affineTransform[x, y, a, b, c, d, e, f] := [a x + c y + e, b x + d y + f] /** Performs an affineTransform with the coefficients specified as a list. */ affineTransform[x, y, coeffs] := {    [a,b,c,d,e,f] = coeffs    return affineTransform[x, y, a, b, c, d, e, f] } /** This performs the inverse of the arbitrary affine transformation above,     converting the point [x',y'] back into the original coordinates [x,y]     using the same affine transform coefficients a,b,c,d,e,f as passed to     the affineTransform function above.     returns:  [x, y] */ inverseAffineTransform[xp, yp, a, b, c, d, e, f] := {    determinant = a d - b c    if determinant == 0       return [undef, undef]    x =   (c (f - yp) - d (e - xp)) / determinant    y =  -(a (f - yp) - b (e - xp)) / determinant        return [x, y] } /* This performs the inverse of the arbitrary affine transformation above,     converting the point [x',y'] back into the original coordinates [x,y]     using the same affine transform coefficients coeffs=[a,b,c,d,e,f] as     passed to the affineTransform function above.     returns:  [x, y] */ inverseAffineTransform[xp, yp, coeffs] := {    [a,b,c,d,e,f] = coeffs    return inverseAffineTransform[xp, yp, a, b, c, d, e, f] } /** This composes a set of transformation coefficients with another set of     transformation coefficients.  The inputs should be an array of 6 elements     containing the coefficients [a,b,c,d,e,f] as described in the documentation     for the affineTransformation function above.     Returns a vector [a,b,c,d,e,f] of the     new coefficients. */ composeTransformation[coeffsA, coeffsB] := {    [a1, b1, c1, d1, e1, f1] = coeffsA    [a2, b2, c2, d2, e2, f2] = coeffsB    return composeTransformation[a1, b1, c1, d1, e1, f1, a2, b2, c2, d2, e2, f2] } /** This composes a set of transformation coefficients with another set of     transformation coefficients.  Returns a vector [a,b,c,d,e,f] of the     new coefficients.  This is just an unrolled, optimized matrix     multiplication, eliminating always-zero terms. */ composeTransformation[a1, b1, c1, d1, e1, f1, a2, b2, c2, d2, e2, f2] := {    return [a1 a2 + b2 c1,        // a            a2 b1 + b2 d1,        // b            a1 c2 + c1 d2,        // c            b1 c2 + d1 d2,        // d            a1 e2 + c1 f2 + e1,   // e            b1 e2 + d1 f2 + f1]   // f } /** Makes a set of transformation coefficients that correspond to a translation     by dx, dy.     This corresponds to the matrix:     x'    a=1  c=0  e=dx     x     y' =  b=0  d=1  f=dy  *  y     1       0    0    1      1     which gives the equation:     x' = x + dx     y' = y + dy */ makeTranslate[dx, dy] := {    return [1,   // a            0,   // b            0,   // c            1,   // d            dx,  // e            dy]  // f } /** Create a scaling transformation that scales the coordinate system     by sx, sy around 0,0.  This corresponds to the matrix:     x'    a=sx  c=0   e=0     x     y' =  b=0   d=sy  f=0  *  y     1       0     0     1     1     which gives the equation:     x' = x sx     y' = y sy */ makeScale[sx, sy] := {    return [sx,  // a            0,   // b            0,   // c            sy,  // d            0,   // e            0]   // f } /** Create a rotating transformation that rotates the coordinate system     around 0,0 by the angle theta.  This corresponds to the matrix:     x'    a=cos[theta]  c=-sin[theta]   e=0     x     y' =  b=sin[theta]  d= cos[theta]   f=0  *  y     1       0              0              1     1     which gives the equations:     x' = x cos[theta] - y sin[theta]     y' = x sin[theta] + y cos[theta] */ makeRotate[theta] := {    c = cos[theta]    s = sin[theta]    return [c,   // a            s,   // b            -s,  // c            c,   // d            0,   // e            0]   // f } /** Create a rotating transformation that rotates the coordinate system     around cx, cy by the angle theta.  This is equivalent to:     translate[cx,cy]     rotate[theta]     translate[-cx, -cy]     This corresponds to the matrix: x'    a=cos[theta]  c=-sin[theta]  e=cx (1 - cos[theta]) + cy sin[theta]     x y' =  b=sin[theta]  d= cos[theta]  f=cy (1 - cos[theta]) - cx sin[theta] *   y  1       0              0             1                                       1     which gives the equations:     x' = (x - cx) cos[theta] - (y - cy) sin[theta] + cx     y' = (x - cx) sin[theta] + (y - cy) cos[theta] + cy */ makeRotateAround[cx, cy, theta] := {    transform = makeTranslate[cx, cy]    transform = composeTransformation[transform, makeRotate[theta]]    transform = composeTransformation[transform, makeTranslate[-cx, -cy]]    return transform } /** Create a scaling transformation that scales the coordinate system     around cx, cy by the scales sx, sy.  This is equivalent to:     translate[cx, cy]     scale[sx, sy]     translate[-cx, -cy]     This corresponds to the matrix:     x'    a=sx   c=0    e=cx (1-sx)        x     y' =  b=0    d=sy   f=cy (1-sy)    *   y      1       0      0      1                1     which gives the equations:     x' = sx (x-cx) + cx     y' = sy (y-cy) + cy    or    x' = sx x + cx (1 - sx)    y' = sy y + cy (1 - sy) */ makeScaleAround[cx, cy, sx, sy] := {    transform = makeTranslate[cx, cy]    transform = composeTransformation[transform, makeScale[sx, sy]]    transform = composeTransformation[transform, makeTranslate[-cx, -cy]]    return transform } /** Make an arbitrary transform with the specified coordinates. */ makeArbitrary[a, b, c, d, e, f] := {    return [a, b, c, d, e, f] } /** Creates the inverse transform of a general transform.  For this to be     valid, the determinant (a d - b c) must be nonzero. */ makeInverse[a, b, c, d, e, f] := {    det = 1 / (a d - b c)    return [ d det,            -b det,            -c det,             a det,            (c f - d e) det,            (b e - a f) det] } /** Creates the inverse transform of a general transform.  For this to be     valid, the determinant (a d - b c) must be nonzero. */ makeInverse[array] := {    [a,b,c,d,e,f] = array    makeInverse[a,b,c,d,e,f] } "Geometry.frink included OK" ```

This is a program written in the programming language Frink.
For more information, view the Frink Documentation or see More Sample Frink Programs.

Alan Eliasen was born 18353 days, 13 hours, 55 minutes ago.