/*
 * topic Bezier patch library
 *
 * This library provides effcient evaluation of
 * Cubic Bezier patches
 */

/**
 * \ingroup rtbezpatch
 * \page rtbezpatch RtBezier Toolkit Overview
 *
 * Overview  \par
 *
 * The toolkit is a group of ten void functions in RenderWare Graphics that 
 * support the Bezier patch code. The functions are useful beyond the 
 * facilities that RenderWare supports, so they are provided and documented in
 * the Toolkit for the developer to use. They use three data types in addition
 * to those documented in the Bezier Patch plugin, RpPatch. \par

 * -# \ref RwV4d is the RenderWare vector that holds four values. Three of them 
 * would normally be x, y and z.
 * -# \ref RtBezierRow is an array of four \ref RwV4d s, one for each control 
 * point of a Bezier spline. Thus any row of control points across a patch can
 * be stored in an \ref RtBezierRow.
 * -# \ref RtBezierMatrix is an array of four \ref RtBezierRow s. So the matrix
 * holds all sixteen control points of a patch. It can hold the coordinates of
 * quad patches or tri patches, and the vectors could represent surface points 
 * instead of control points. Alternatively it can represent normal vectors,
 * tangents and "difference" values and all these are used in the Bezier 
 * Toolkit functions. \par 

 *
 * A Bezier curve can be measured in "u". At the first control point u=0.
 * At the last, u=1 and at the other two control points, u=1/3 and u=2/3. This
 * is applied to patches, so that one edge is measured by u. The transverse
 * edge is also a Bezier curved measured in "v", with control points at v=0,
 * v=1/3, v=2/3 and v=1. This notation is used in some functions.
 *
 * The Tool kit functions fall into four groups. \par
 *
 * -# \ref RtBezierQuadFromTriangle() is given only the data for a tri patch
 * and it writes an \ref RtBezierMatrix containing the control points of a complete
 * Bezier quad patch. This is a utility function to serve several other 
 * functions that take only quad patches as parameters.
 * -# \ref RtBezierQuadMatrixControlFit(),
 * \ref RtBezierTriMatrixControlFit() and 
 * \ref RtBezierQuadMatrixSample() convert from control points to surface 
 * points. The first two calculate the surface point that corresponds to each
 * control point of tri and of quad patches. The third will calculate any point
 * on the patch's surface.
 * -# \ref RtBezierQuadMatrixWeightSetup() pre-calculates many of the 
 * matrix-multiplications needed for patches. 
 * \ref RtBezierQuadMatrixDifferenceSetup() calculates values for "forward 
 * differencing" to speed up the calculation of patch coordinates. They are 
 * used with \ref RtBezierQuadDifferenceStepU() and 
 * \ref RtBezierQuadDifferenceStepV(), which calculate values for rows of
 * vertices and for vertices over a whole patch.
 * -# \ref RtBezierQuadTangent() and \ref RtBezierQuadMatrixGetNormals() 
 * calculate the directions of normal vectors and of tangents at each vertex of
 * a patch.
 *
 */

#include <stdarg.h>
#include <stddef.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>

#include <rwcore.h>
#include <rpdbgerr.h>

#include "rtbezpat.h"

static const char   rcsid[] __RWUNUSED__ =
    "@@(#)$Id: rtbezpat.c,v 1.22 2001/10/05 11:42:34 rburg Exp $";

#define  RtBezierQuadFromTriangleMacro(P, T)                             \
MACRO_START                                                              \
{                                                                        \
    const RtBezierV4d  *const T0 = &(T)[0][0];                           \
    const RtBezierV4d  *const T1 = &(T)[1][0];                           \
    const RtBezierV4d  *const T2 = &(T)[2][0];                           \
    const RtBezierV4d  *const T3 = &(T)[3][0];                           \
    RtBezierV4d        *const P0 = &(P)[0][0];                           \
    RtBezierV4d        *const P1 = &(P)[1][0];                           \
    RtBezierV4d        *const P2 = &(P)[2][0];                           \
    RtBezierV4d        *const P3 = &(P)[3][0];                           \
                                                                         \
                                                                         \
    /* Row 0 */                                                          \
    P0[0].x = (T0[0].x);                                                 \
    P0[0].y = (T0[0].y);                                                 \
    P0[0].z = (T0[0].z);                                                 \
    P0[1].x = (T0[1].x);                                                 \
    P0[1].y = (T0[1].y);                                                 \
    P0[1].z = (T0[1].z);                                                 \
    P0[2].x = (T0[2].x);                                                 \
    P0[2].y = (T0[2].y);                                                 \
    P0[2].z = (T0[2].z);                                                 \
    P0[3].x = (T0[3].x);                                                 \
    P0[3].y = (T0[3].y);                                                 \
    P0[3].z = (T0[3].z);                                                 \
                                                                         \
    /* Row 1 */                                                          \
    P1[0].x = (T1[0].x);                                                 \
    P1[0].y = (T1[0].y);                                                 \
    P1[0].z = (T1[0].z);                                                 \
    P1[1].x = ((2 * T1[1].x + T1[0].x - T0[0].x + T0[1].x) / 3);         \
    P1[1].y = ((2 * T1[1].y + T1[0].y - T0[0].y + T0[1].y) / 3);         \
    P1[1].z = ((2 * T1[1].z + T1[0].z - T0[0].z + T0[1].z) / 3);         \
    P1[2].x = ((T1[2].x + 2 * T1[1].x + 2 * T0[2].x - 2 * T0[1].x) / 3); \
    P1[2].y = ((T1[2].y + 2 * T1[1].y + 2 * T0[2].y - 2 * T0[1].y) / 3); \
    P1[2].z = ((T1[2].z + 2 * T1[1].z + 2 * T0[2].z - 2 * T0[1].z) / 3); \
    P1[3].x = (T1[2].x - T0[2].x + T0[3].x);                             \
    P1[3].y = (T1[2].y - T0[2].y + T0[3].y);                             \
    P1[3].z = (T1[2].z - T0[2].z + T0[3].z);                             \
                                                                         \
    /* Row 2 */                                                          \
    P2[0].x = (T2[0].x);                                                 \
    P2[0].y = (T2[0].y);                                                 \
    P2[0].z = (T2[0].z);                                                 \
    P2[1].x = ((T2[1].x + 2 * T2[0].x + 2 * T1[1].x - 2 * T1[0].x) / 3); \
    P2[1].y = ((T2[1].y + 2 * T2[0].y + 2 * T1[1].y - 2 * T1[0].y) / 3); \
    P2[1].z = ((T2[1].z + 2 * T2[0].z + 2 * T1[1].z - 2 * T1[0].z) / 3); \
    P2[2].x = ((2 * T2[1].x + T2[0].x + 2 * T1[2].x - 2 * T1[0].x +      \
                T0[0].x + T0[2].x - 2 * T0[1].x) / 3);                   \
    P2[2].y = ((2 * T2[1].y + T2[0].y + 2 * T1[2].y - 2 * T1[0].y +      \
                T0[0].y + T0[2].y - 2 * T0[1].y) / 3);                   \
    P2[2].z = ((2 * T2[1].z + T2[0].z + 2 * T1[2].z - 2 * T1[0].z +      \
                T0[0].z + T0[2].z - 2 * T0[1].z) / 3);                   \
    P2[3].x = (T0[1].x - 2 * T1[1].x + T2[1].x + 2 * T1[2].x -           \
               2 * T0[2].x + T0[3].x);                                   \
    P2[3].y = (T0[1].y - 2 * T1[1].y + T2[1].y + 2 * T1[2].y -           \
               2 * T0[2].y + T0[3].y);                                   \
    P2[3].z = (T0[1].z - 2 * T1[1].z + T2[1].z + 2 * T1[2].z -           \
               2 * T0[2].z + T0[3].z);                                   \
                                                                         \
    /* Row 3 */                                                          \
    P3[0].x = (T3[0].x);                                                 \
    P3[0].y = (T3[0].y);                                                 \
    P3[0].z = (T3[0].z);                                                 \
    P3[1].x = (T3[0].x - T2[0].x + T2[1].x);                             \
    P3[1].y = (T3[0].y - T2[0].y + T2[1].y);                             \
    P3[1].z = (T3[0].z - T2[0].z + T2[1].z);                             \
    P3[2].x = (T3[0].x - 2 * T2[0].x + 2 * T2[1].x + T1[2].x -           \
               2 * T1[1].x + T1[0].x);                                   \
    P3[2].y = (T3[0].y - 2 * T2[0].y + 2 * T2[1].y + T1[2].y -           \
               2 * T1[1].y + T1[0].y);                                   \
    P3[2].z = (T3[0].z - 2 * T2[0].z + 2 * T2[1].z + T1[2].z -           \
               2 * T1[1].z + T1[0].z);                                   \
    P3[3].x = (T3[0].x + 3 * T2[1].x - 3 * T2[0].x + 3 * T1[2].x -       \
               6 * T1[1].x + 3 * T1[0].x + T0[3].x - 3 * T0[2].x +       \
               3 * T0[1].x - T0[0].x);                                   \
    P3[3].y = (T3[0].y + 3 * T2[1].y - 3 * T2[0].y + 3 * T1[2].y -       \
               6 * T1[1].y + 3 * T1[0].y + T0[3].y - 3 * T0[2].y +       \
               3 * T0[1].y - T0[0].y);                                   \
    P3[3].z = (T3[0].z + 3 * T2[1].z - 3 * T2[0].z + 3 * T1[2].z -       \
               6 * T1[1].z + 3 * T1[0].z + T0[3].z - 3 * T0[2].z +       \
               3 * T0[1].z - T0[0].z);                                   \
}                                                                        \
MACRO_STOP

/**
 * \ingroup rtbezpatch
 * \ref RtBezierQuadFromTriangle takes a matrix of 10 control points for 
 * a tri patch, and returns a matrix of the 16 control points for a 
 * quad patch. The side of the tri patch defined by control points 9-8-6-3
 * becomes quad patch's diagonal 12-9-6-3.
 *
 * \param P   Output control points of Bezier quadrilateral
 * \param T   Input  control points of Bezier triangle
 *
 */
void
RtBezierQuadFromTriangle(RtBezierMatrix P, RtBezierMatrix T)
{
    RWAPIFUNCTION(RWSTRING("RtBezierQuadFromTriangle"));

    RtBezierQuadFromTriangleMacro(P, T);

    RWRETURNVOID();
}
