#include "allheaders.h"
/*-------------------------------------------------------------*
* Composable coordinate transforms *
*-------------------------------------------------------------*/
/*!
* \brief createMatrix2dTranslate()
*
* \param[in] transx x component of translation wrt. the origin
* \param[in] transy y component of translation wrt. the origin
* \return 3x3 transform matrix, or NULL on error
*
*
* Notes:
* (1) The translation is equivalent to:
* v' = Av
* where v and v' are 1x3 column vectors in the form
* v = [x, y, 1]^ ^ denotes transpose
* and the affine translation matrix is
* A = [ 1 0 tx
* 0 1 ty
* 0 0 1 ]
*
* (2) We consider translation as with respect to a fixed origin.
* In a clipping operation, the origin moves and the points
* are fixed, and you use (-tx, -ty) where (tx, ty) is the
* translation vector of the origin.
*
*/
l_float32 *
createMatrix2dTranslate(l_float32 transx,
l_float32 transy)
{
l_float32 *mat;
mat = (l_float32 *)LEPT_CALLOC(9, sizeof(l_float32));
mat[0] = mat[4] = mat[8] = 1;
mat[2] = transx;
mat[5] = transy;
return mat;
}
/*!
* \brief createMatrix2dScale()
*
* \param[in] scalex horizontal scale factor
* \param[in] scaley vertical scale factor
* \return 3x3 transform matrix, or NULL on error
*
*
* Notes:
* (1) The scaling is equivalent to:
* v' = Av
* where v and v' are 1x3 column vectors in the form
* v = [x, y, 1]^ ^ denotes transpose
* and the affine scaling matrix is
* A = [ sx 0 0
* 0 sy 0
* 0 0 1 ]
*
* (2) We consider scaling as with respect to a fixed origin.
* In other words, the origin is the only point that doesn't
* move in the scaling transform.
*
*/
l_float32 *
createMatrix2dScale(l_float32 scalex,
l_float32 scaley)
{
l_float32 *mat;
mat = (l_float32 *)LEPT_CALLOC(9, sizeof(l_float32));
mat[0] = scalex;
mat[4] = scaley;
mat[8] = 1;
return mat;
}
/*!
* \brief createMatrix2dRotate()
*
* \param[in] xc, yc location of center of rotation
* \param[in] angle rotation in radians; clockwise is positive
* \return 3x3 transform matrix, or NULL on error
*
*
* Notes:
* (1) The rotation is equivalent to:
* v' = Av
* where v and v' are 1x3 column vectors in the form
* v = [x, y, 1]^ ^ denotes transpose
* and the affine rotation matrix is
* A = [ cosa -sina xc*1-cosa + yc*sina
* sina cosa yc*1-cosa - xc*sina
* 0 0 1 ]
*
* If the rotation is about the origin, xc, yc) = (0, 0 and
* this simplifies to
* A = [ cosa -sina 0
* sina cosa 0
* 0 0 1 ]
*
* These relations follow from the following equations, which
* you can convince yourself are correct as follows. Draw a
* circle centered on xc,yc) and passing through (x,y), with
* (x',y') on the arc at an angle 'a' clockwise from (x,y).
* [ Hint: cosa + b = cosa * cosb - sina * sinb
* sina + b = sina * cosb + cosa * sinb ]
*
* x' - xc = x - xc) * cosa - (y - yc * sina
* y' - yc = x - xc) * sina + (y - yc * cosa
*
*/
l_float32 *
createMatrix2dRotate(l_float32 xc,
l_float32 yc,
l_float32 angle)
{
l_float32 sina, cosa;
l_float32 *mat;
mat = (l_float32 *)LEPT_CALLOC(9, sizeof(l_float32));
sina = sin(angle);
cosa = cos(angle);
mat[0] = mat[4] = cosa;
mat[1] = -sina;
mat[2] = xc * (1.0 - cosa) + yc * sina;
mat[3] = sina;
mat[5] = yc * (1.0 - cosa) - xc * sina;
mat[8] = 1;
return mat;
}
/*-------------------------------------------------------------*
* Special coordinate transforms on pta *
*-------------------------------------------------------------*/
/*!
* \brief ptaTranslate()
*
* \param[in] ptas for initial points
* \param[in] transx x component of translation wrt. the origin
* \param[in] transy y component of translation wrt. the origin
* \return ptad translated points, or NULL on error
*
*
* Notes:
* (1) See createMatrix2dTranslate() for details of transform.
*
*/
PTA *
ptaTranslate(PTA *ptas,
l_float32 transx,
l_float32 transy)
{
l_int32 i, npts;
l_float32 x, y;
PTA *ptad;
PROCNAME("ptaTranslate");
if (!ptas)
return (PTA *)ERROR_PTR("ptas not defined", procName, NULL);
npts = ptaGetCount(ptas);
if ((ptad = ptaCreate(npts)) == NULL)
return (PTA *)ERROR_PTR("ptad not made", procName, NULL);
for (i = 0; i < npts; i++) {
ptaGetPt(ptas, i, &x, &y);
ptaAddPt(ptad, x + transx, y + transy);
}
return ptad;
}
/*!
* \brief ptaScale()
*
* \param[in] ptas for initial points
* \param[in] scalex horizontal scale factor
* \param[in] scaley vertical scale factor
* \return 0 if OK; 1 on error
*
*
* Notes:
* (1) See createMatrix2dScale() for details of transform.
*
*/
PTA *
ptaScale(PTA *ptas,
l_float32 scalex,
l_float32 scaley)
{
l_int32 i, npts;
l_float32 x, y;
PTA *ptad;
PROCNAME("ptaScale");
if (!ptas)
return (PTA *)ERROR_PTR("ptas not defined", procName, NULL);
npts = ptaGetCount(ptas);
if ((ptad = ptaCreate(npts)) == NULL)
return (PTA *)ERROR_PTR("ptad not made", procName, NULL);
for (i = 0; i < npts; i++) {
ptaGetPt(ptas, i, &x, &y);
ptaAddPt(ptad, scalex * x, scaley * y);
}
return ptad;
}
/*!
* \brief ptaRotate()
*
* \param[in] ptas for initial points
* \param[in] xc, yc location of center of rotation
* \param[in] angle rotation in radians; clockwise is positive
* \return 0 if OK; 1 on error
*
*
* Notes;
* (1) See createMatrix2dScale() for details of transform.
* (2) This transform can be thought of as composed of the
* sum of two parts:
* a) an (x,y)-dependent rotation about the origin:
* xr = x * cosa - y * sina
* yr = x * sina + y * cosa
* b) an (x,y)-independent translation that depends on the
* rotation center and the angle:
* xt = xc - xc * cosa + yc * sina
* yt = yc - xc * sina - yc * cosa
* The translation part (xt,yt) is equal to the difference
* between the center (xc,yc) and the location of the
* center after it is rotated about the origin.
*
*/
PTA *
ptaRotate(PTA *ptas,
l_float32 xc,
l_float32 yc,
l_float32 angle)
{
l_int32 i, npts;
l_float32 x, y, xp, yp, sina, cosa;
PTA *ptad;
PROCNAME("ptaRotate");
if (!ptas)
return (PTA *)ERROR_PTR("ptas not defined", procName, NULL);
npts = ptaGetCount(ptas);
if ((ptad = ptaCreate(npts)) == NULL)
return (PTA *)ERROR_PTR("ptad not made", procName, NULL);
sina = sin(angle);
cosa = cos(angle);
for (i = 0; i < npts; i++) {
ptaGetPt(ptas, i, &x, &y);
xp = xc + (x - xc) * cosa - (y - yc) * sina;
yp = yc + (x - xc) * sina + (y - yc) * cosa;
ptaAddPt(ptad, xp, yp);
}
return ptad;
}
/*-------------------------------------------------------------*
* Special coordinate transforms on boxa *
*-------------------------------------------------------------*/
/*!
* \brief boxaTranslate()
*
* \param[in] boxas
* \param[in] transx x component of translation wrt. the origin
* \param[in] transy y component of translation wrt. the origin
* \return boxad translated boxas, or NULL on error
*
* Notes:
* (1) See createMatrix2dTranslate() for details of transform.
*/
BOXA *
boxaTranslate(BOXA *boxas,
l_float32 transx,
l_float32 transy)
{
PTA *ptas, *ptad;
BOXA *boxad;
PROCNAME("boxaTranslate");
if (!boxas)
return (BOXA *)ERROR_PTR("boxas not defined", procName, NULL);
ptas = boxaConvertToPta(boxas, 4);
ptad = ptaTranslate(ptas, transx, transy);
boxad = ptaConvertToBoxa(ptad, 4);
ptaDestroy(&ptas);
ptaDestroy(&ptad);
return boxad;
}
/*!
* \brief boxaScale()
*
* \param[in] boxas
* \param[in] scalex horizontal scale factor
* \param[in] scaley vertical scale factor
* \return boxad scaled boxas, or NULL on error
*
* Notes:
* (1) See createMatrix2dScale() for details of transform.
*/
BOXA *
boxaScale(BOXA *boxas,
l_float32 scalex,
l_float32 scaley)
{
PTA *ptas, *ptad;
BOXA *boxad;
PROCNAME("boxaScale");
if (!boxas)
return (BOXA *)ERROR_PTR("boxas not defined", procName, NULL);
ptas = boxaConvertToPta(boxas, 4);
ptad = ptaScale(ptas, scalex, scaley);
boxad = ptaConvertToBoxa(ptad, 4);
ptaDestroy(&ptas);
ptaDestroy(&ptad);
return boxad;
}
/*!
* \brief boxaRotate()
*
* \param[in] boxas
* \param[in] xc, yc location of center of rotation
* \param[in] angle rotation in radians; clockwise is positive
* \return boxad scaled boxas, or NULL on error
*
* Notes:
* (1) See createMatrix2dRotate() for details of transform.
*/
BOXA *
boxaRotate(BOXA *boxas,
l_float32 xc,
l_float32 yc,
l_float32 angle)
{
PTA *ptas, *ptad;
BOXA *boxad;
PROCNAME("boxaRotate");
if (!boxas)
return (BOXA *)ERROR_PTR("boxas not defined", procName, NULL);
ptas = boxaConvertToPta(boxas, 4);
ptad = ptaRotate(ptas, xc, yc, angle);
boxad = ptaConvertToBoxa(ptad, 4);
ptaDestroy(&ptas);
ptaDestroy(&ptad);
return boxad;
}
/*-------------------------------------------------------------*
* General affine coordinate transform *
*-------------------------------------------------------------*/
/*!
* \brief ptaAffineTransform()
*
* \param[in] ptas for initial points
* \param[in] mat 3x3 transform matrix; canonical form
* \return ptad transformed points, or NULL on error
*/
PTA *
ptaAffineTransform(PTA *ptas,
l_float32 *mat)
{
l_int32 i, npts;
l_float32 vecs[3], vecd[3];
PTA *ptad;
PROCNAME("ptaAffineTransform");
if (!ptas)
return (PTA *)ERROR_PTR("ptas not defined", procName, NULL);
if (!mat)
return (PTA *)ERROR_PTR("transform not defined", procName, NULL);
vecs[2] = 1;
npts = ptaGetCount(ptas);
if ((ptad = ptaCreate(npts)) == NULL)
return (PTA *)ERROR_PTR("ptad not made", procName, NULL);
for (i = 0; i < npts; i++) {
ptaGetPt(ptas, i, &vecs[0], &vecs[1]);
l_productMatVec(mat, vecs, vecd, 3);
ptaAddPt(ptad, vecd[0], vecd[1]);
}
return ptad;
}
/*!
* \brief boxaAffineTransform()
*
* \param[in] boxas
* \param[in] mat 3x3 transform matrix; canonical form
* \return boxad transformed boxas, or NULL on error
*/
BOXA *
boxaAffineTransform(BOXA *boxas,
l_float32 *mat)
{
PTA *ptas, *ptad;
BOXA *boxad;
PROCNAME("boxaAffineTransform");
if (!boxas)
return (BOXA *)ERROR_PTR("boxas not defined", procName, NULL);
if (!mat)
return (BOXA *)ERROR_PTR("transform not defined", procName, NULL);
ptas = boxaConvertToPta(boxas, 4);
ptad = ptaAffineTransform(ptas, mat);
boxad = ptaConvertToBoxa(ptad, 4);
ptaDestroy(&ptas);
ptaDestroy(&ptad);
return boxad;
}
/*-------------------------------------------------------------*
* Matrix operations *
*-------------------------------------------------------------*/
/*!
* \brief l_productMatVec()
*
* \param[in] mat square matrix, as a 1-dimensional %size^2 array
* \param[in] vecs input column vector of length %size
* \param[in] vecd result column vector
* \param[in] size matrix is %size x %size; vectors are length %size
* \return 0 if OK, 1 on error
*/
l_ok
l_productMatVec(l_float32 *mat,
l_float32 *vecs,
l_float32 *vecd,
l_int32 size)
{
l_int32 i, j;
PROCNAME("l_productMatVec");
if (!mat)
return ERROR_INT("matrix not defined", procName, 1);
if (!vecs)
return ERROR_INT("input vector not defined", procName, 1);
if (!vecd)
return ERROR_INT("result vector not defined", procName, 1);
for (i = 0; i < size; i++) {
vecd[i] = 0;
for (j = 0; j < size; j++) {
vecd[i] += mat[size * i + j] * vecs[j];
}
}
return 0;
}
/*!
* \brief l_productMat2()
*
* \param[in] mat1 square matrix, as a 1-dimensional size^2 array
* \param[in] mat2 square matrix, as a 1-dimensional size^2 array
* \param[in] matd square matrix; product stored here
* \param[in] size of matrices
* \return 0 if OK, 1 on error
*/
l_ok
l_productMat2(l_float32 *mat1,
l_float32 *mat2,
l_float32 *matd,
l_int32 size)
{
l_int32 i, j, k, index;
PROCNAME("l_productMat2");
if (!mat1)
return ERROR_INT("matrix 1 not defined", procName, 1);
if (!mat2)
return ERROR_INT("matrix 2 not defined", procName, 1);
if (!matd)
return ERROR_INT("result matrix not defined", procName, 1);
for (i = 0; i < size; i++) {
for (j = 0; j < size; j++) {
index = size * i + j;
matd[index] = 0;
for (k = 0; k < size; k++)
matd[index] += mat1[size * i + k] * mat2[size * k + j];
}
}
return 0;
}
/*!
* \brief l_productMat3()
*
* \param[in] mat1 square matrix, as a 1-dimensional size^2 array
* \param[in] mat2 square matrix, as a 1-dimensional size^2 array
* \param[in] mat3 square matrix, as a 1-dimensional size^2 array
* \param[in] matd square matrix; product stored here
* \param[in] size of matrices
* \return 0 if OK, 1 on error
*/
l_ok
l_productMat3(l_float32 *mat1,
l_float32 *mat2,
l_float32 *mat3,
l_float32 *matd,
l_int32 size)
{
l_float32 *matt;
PROCNAME("l_productMat3");
if (!mat1)
return ERROR_INT("matrix 1 not defined", procName, 1);
if (!mat2)
return ERROR_INT("matrix 2 not defined", procName, 1);
if (!mat3)
return ERROR_INT("matrix 3 not defined", procName, 1);
if (!matd)
return ERROR_INT("result matrix not defined", procName, 1);
if ((matt = (l_float32 *)LEPT_CALLOC((size_t)size * size,
sizeof(l_float32))) == NULL)
return ERROR_INT("matt not made", procName, 1);
l_productMat2(mat1, mat2, matt, size);
l_productMat2(matt, mat3, matd, size);
LEPT_FREE(matt);
return 0;
}
/*!
* \brief l_productMat4()
*
* \param[in] mat1 square matrix, as a 1-dimensional size^2 array
* \param[in] mat2 square matrix, as a 1-dimensional size^2 array
* \param[in] mat3 square matrix, as a 1-dimensional size^2 array
* \param[in] mat4 square matrix, as a 1-dimensional size^2 array
* \param[in] matd square matrix; product stored here
* \param[in] size of matrices
* \return 0 if OK, 1 on error
*/
l_ok
l_productMat4(l_float32 *mat1,
l_float32 *mat2,
l_float32 *mat3,
l_float32 *mat4,
l_float32 *matd,
l_int32 size)
{
l_float32 *matt;
PROCNAME("l_productMat4");
if (!mat1)
return ERROR_INT("matrix 1 not defined", procName, 1);
if (!mat2)
return ERROR_INT("matrix 2 not defined", procName, 1);
if (!mat3)
return ERROR_INT("matrix 3 not defined", procName, 1);
if (!matd)
return ERROR_INT("result matrix not defined", procName, 1);
if ((matt = (l_float32 *)LEPT_CALLOC((size_t)size * size,
sizeof(l_float32))) == NULL)
return ERROR_INT("matt not made", procName, 1);
l_productMat3(mat1, mat2, mat3, matt, size);
l_productMat2(matt, mat4, matd, size);
LEPT_FREE(matt);
return 0;
}