scilab/modules/graphics/src/c/math_graphics.c

   1 /*
   2  * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
   3  * Copyright (C) 1998-2001 - ENPC - Jean-Philippe Chancelier
   4  * Copyright (C) 2006 - INRIA - Jean-Baptiste Silvy
   5  *
   6  * This file must be used under the terms of the CeCILL.
   7  * This source file is licensed as described in the file COPYING, which
   8  * you should have received as part of this distribution.  The terms
   9  * are also available at
  10  * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
  11  *
  12  */
  13
  14 /*------------------------------------------------------------------------
  15  *    Graphic library
  16  --------------------------------------------------------------------------*/
  17
  18 #include "math_graphics.h"
  19
  20 /*
  21  * we use spConfig.h for machine constants
  22  * XXX : spConfig should be merged and unified
  23  *       with other machine constant scilab code
  24  */
  25
  26 #define spINSIDE_SPARSE
  27 #include "../../sparse/includes/spConfig.h"
  28
  29 double Mini(const double vect[], int n)
  30 {
  31   int i;
  32   double vmin;
  33   vmin = LARGEST_REAL;
  34   for (i = 0 ; i < n ; i++)
  35     /*    if ( isinf(vect[i])== 0 && isnan(vect[i])==0 && vect[i] < vmin)  */
  36     if ( finite(vect[i])== 1 && vect[i] < vmin)
  37       vmin=vect[i];
  38   return(vmin);
  39 }
  40
  41 double Maxi(const double vect[], int n)
  42 {
  43   int i;
  44   double maxi;
  45   maxi= - LARGEST_REAL;
  46   for (i =0 ; i < n ; i++)
  47     /* if ( isinf(vect[i])== 0 && isnan(vect[i])==0 && vect[i] > maxi) */
  48     if ( finite(vect[i])== 1 && vect[i] > maxi)
  49       maxi=vect[i];
  50   return(maxi);
  51 }
  52
  53 /*----------------------------------------------------------------------------*/
  54
  55 /* perform the rotation of point from to point dest  */
  56 void rotate2D( double from[2], double center[2], double angle, double dest[2] )
  57 {
  58   rotate2Dim( from, center, cos( angle ), sin( angle ), dest ) ;
  59 }
  60
  61 /*----------------------------------------------------------------------------*/
  62 /* perform the rotation of point from to point to. */
  63 /* the angle is directly given with its sine and cosine for speed */
  64 void rotate2Dim( double from[2]   ,
  65                  double center[2] ,
  66                  double cosAngle  ,
  67                  double sinAngle  ,
  68                  double dest[2]    )
  69 {
  70   double diff[2] ;
  71
  72   /* put the center to (0,0) */
  73   diff[0] = from[0] - center[0] ;
  74   diff[1] = from[1] - center[1] ;
  75
  76   /* turn and translate back */
  77   dest[0] = diff[0] * cosAngle - diff[1] * sinAngle + center[0] ;
  78   dest[1] = diff[0] * sinAngle + diff[1] * cosAngle + center[1] ;
  79 }
  80 /*----------------------------------------------------------------------------*/
  81 /* perform the translation of point from to point to with vector trans */
  82 void translate2D( double from[2], double trans[2], double dest[2] )
  83 {
  84   dest[0] = from[0] + trans[0] ;
  85   dest[1] = from[1] + trans[1] ;
  86 }
  87 /*----------------------------------------------------------------------------*/
  88 void iTranslate2D( int from[2], int trans[2], int dest[2] )
  89 {
  90   dest[0] = from[0] + trans[0] ;
  91   dest[1] = from[1] + trans[1] ;
  92 }
  93 /*----------------------------------------------------------------------------*/
  94 void vectSubstract2D(const double vect1[2], const double vect2[], double res[2])
  95 {
  96   res[0] = vect1[0] - vect2[0];
  97   res[1] = vect1[1] - vect2[1];
  98 }
  99 /*----------------------------------------------------------------------------*/
 100 void vectAdd2D(const double v1[2], const double v2[2], double res[2])
 101 {
 102   res[0] = v1[0] + v2[0];
 103   res[1] = v1[1] + v2[1];
 104 }
 105 /*----------------------------------------------------------------------------*/
 106 void scalarMult2D(const double v[2], const double scalar, double res[2])
 107 {
 108   res[0] = scalar * v[0];
 109   res[1] = scalar * v[1];
 110 }
 111 /*----------------------------------------------------------------------------*/
 112 void normalize2d( double vect[2] )
 113 {
 114   double norm = NORM_2D(vect) ;
 115   vect[0] /= norm ;
 116   vect[1] /= norm ;
 117 }
 118 /*----------------------------------------------------------------------------*/
 119 void iNormalize2d( int vect[2] )
 120 {
 121   double norm = NORM_2D(vect) ;
 122   vect[0] = round( vect[0] / norm ) ;
 123   vect[1] = round( vect[1] / norm ) ;
 124 }
 125 /*----------------------------------------------------------------------------*/
 126 BOOL isPointInTriangle(const double point[2], const double a[2],
 127                        const double b[2], const double c[2])
 128 {
 129   return (   areOnSameSideOfLine(point, a, b, c)
 130           && areOnSameSideOfLine(point, b, a, c)
 131           && areOnSameSideOfLine(point, c, a, b));
 132 }
 133 /*----------------------------------------------------------------------------*/
 134 BOOL areOnSameSideOfLine(const double p1[2], const double p2[2],
 135                          const double a[2], const double b[2])
 136 {
 137   // point are on the same if and only if (AB^AP1).(AB^AP2) >= 0
 138   double ab[3];
 139   double ap1[3];
 140   double ap2[3];
 141   double cp1[3];
 142   double cp2[3];
 143
 144   ab[0] = b[0] - a[0];
 145   ab[1] = b[1] - a[1];
 146   ab[2] = 0.0;
 147
 148   ap1[0] = p1[0] - a[0];
 149   ap1[1] = p1[1] - a[1];
 150   ap1[2] = 0.0;
 151
 152   ap2[0] = p2[0] - a[0];
 153   ap2[1] = p2[1] - a[1];
 154   ap2[2] = 0.0;
 155
 156   crossProduct(ab, ap1, cp1);
 157   crossProduct(ab, ap2, cp2);
 158
 159   return (DOT_PROD_3D(cp1, cp2) >= 0.0);
 160
 161 }
 162 /*----------------------------------------------------------------------------*/
 163 void crossProduct( const double v1[3], const double v2[3], double res[3] )
 164 {
 165   /* save data to be able to use v1 o v2 as res */
 166   double v10 = v1[0];
 167   double v20 = v2[0];
 168   double v11 = v1[1];
 169   double v21 = v2[1];
 170   res[0] = v11 * v2[2] - v1[2] * v21 ;
 171   res[1] = v1[2] * v20 - v10 * v2[2] ;
 172   res[2] = v10 * v21 - v11 * v20 ;
 173 }
 174 /*----------------------------------------------------------------------------*/
 175 void vectSubstract3D(const double v1[3] ,const double v2[3], double res[3])
 176 {
 177   res[0] = v1[0] - v2[0];
 178   res[1] = v1[1] - v2[1];
 179   res[2] = v1[2] - v2[2];
 180 }
 181 /*----------------------------------------------------------------------------*/
 182 void vectAdd3D(const double v1[3], const double v2[3], double res[3])
 183 {
 184   res[0] = v1[0] + v2[0];
 185   res[1] = v1[1] + v2[1];
 186   res[2] = v1[2] + v2[2];
 187 }
 188 /*----------------------------------------------------------------------------*/
 189 void scalarMult3D(const double v[3], double scalar, double res[3])
 190 {
 191   res[0] = scalar * v[0];
 192   res[1] = scalar * v[1];
 193   res[2] = scalar * v[2];
 194 }
 195 /*----------------------------------------------------------------------------*/
 196 void normalize3D( double vect[3] )
 197 {
 198   double norm = NORM_3D(vect) ;
 199   vect[0] /= norm ;
 200   vect[1] /= norm ;
 201   vect[2] /= norm ;
 202 }
 203 /*----------------------------------------------------------------------------*/
 204 void setToIdentity(double mat4D[4][4])
 205 {
 206   int i;
 207   int j;
 208   for (i = 0; i < 4; i++)
 209   {
 210     for (j = 0; j < 4; j++)
 211     {
 212       mat4D[i][j] = 0.0;
 213     }
 214     mat4D[i][i] = 1.0;
 215   }
 216 }
 217 /*----------------------------------------------------------------------------*/
 218 void mat4DMult(const double mat4D[4][4], const double vect3D[3], double res[3])
 219 {
 220   int i;
 221   double res4D[4];
 222   // w coordinate of the vector is supposed to be 1;
 223   for (i = 0; i < 4; i++) {
 224     res4D[i]  =  vect3D[0] * mat4D[i][0] + vect3D[1] * mat4D[i][1]
 225                + vect3D[2] * mat4D[i][2] + mat4D[i][3];
 226   }
 227   res[0] = res4D[0] / res4D[3];
 228   res[1] = res4D[1] / res4D[3];
 229   res[2] = res4D[2] / res4D[3];
 230 }
 231 /*----------------------------------------------------------------------------*/