main/editor/fixseg.c

   1 /*
   2 THE COMPUTER CODE CONTAINED HEREIN IS THE SOLE PROPERTY OF PARALLAX
   3 SOFTWARE CORPORATION ("PARALLAX").  PARALLAX, IN DISTRIBUTING THE CODE TO
   4 END-USERS, AND SUBJECT TO ALL OF THE TERMS AND CONDITIONS HEREIN, GRANTS A
   5 ROYALTY-FREE, PERPETUAL LICENSE TO SUCH END-USERS FOR USE BY SUCH END-USERS
   6 IN USING, DISPLAYING,  AND CREATING DERIVATIVE WORKS THEREOF, SO LONG AS
   7 SUCH USE, DISPLAY OR CREATION IS FOR NON-COMMERCIAL, ROYALTY OR REVENUE
   8 FREE PURPOSES.  IN NO EVENT SHALL THE END-USER USE THE COMPUTER CODE
   9 CONTAINED HEREIN FOR REVENUE-BEARING PURPOSES.  THE END-USER UNDERSTANDS
  10 AND AGREES TO THE TERMS HEREIN AND ACCEPTS THE SAME BY USE OF THIS FILE.
  11 COPYRIGHT 1993-1998 PARALLAX SOFTWARE CORPORATION.  ALL RIGHTS RESERVED.
  12 */
  13
  14 /*
  15  *
  16  * Functions to make faces planar and probably other things.
  17  *
  18  */
  19
  20 #ifdef HAVE_CONFIG_H
  21 #include "conf.h"
  22 #endif
  23
  24 #include <stdio.h>
  25 #include <stdlib.h>
  26 #include <math.h>
  27 #include <string.h>
  28
  29 #include "mono.h"
  30 #include "key.h"
  31 #include "gr.h"
  32 #include "inferno.h"
  33 #include "editor.h"
  34 #include "dxxerror.h"
  35
  36
  37 #define SWAP(a,b) {temp = (a); (a) = (b); (b) = temp;}
  38
  39 //      -----------------------------------------------------------------------------------------------------------------
  40 //      Gauss-Jordan elimination solution of a system of linear equations.
  41 //      a[1..n][1..n] is the input matrix.  b[1..n][1..m] is input containing the m right-hand side vectors.
  42 //      On output, a is replaced by its matrix inverse and b is replaced by the corresponding set of solution vectors.
  43 void gaussj(fix **a, int n, fix **b, int m)
  44 {
  45         int     indxc[4], indxr[4], ipiv[4];
  46         int     i, icol=0, irow=0, j, k, l, ll;
  47         fix     big, dum, pivinv, temp;
  48
  49         if (n > 4) {
  50                 mprintf((0,"Error -- array too large in gaussj.\n"));
  51                 Int3();
  52         }
  53
  54         for (j=1; j<=n; j++)
  55                 ipiv[j] = 0;
  56
  57         for (i=1; i<=n; i++) {
  58                 big = 0;
  59                 for (j=1; j<=n; j++)
  60                         if (ipiv[j] != 1)
  61                                 for (k=1; k<=n; k++) {
  62                                         if (ipiv[k] == 0) {
  63                                                 if (abs(a[j][k]) >= big) {
  64                                                         big = abs(a[j][k]);
  65                                                         irow = j;
  66                                                         icol = k;
  67                                                 }
  68                                         } else if (ipiv[k] > 1) {
  69                                                 mprintf((0,"Error: Singular matrix-1\n"));
  70                                                 Int3();
  71                                         }
  72                                 }
  73
  74                 ++(ipiv[icol]);
  75
  76                 // We now have the pivot element, so we interchange rows, if needed, to put the pivot
  77                 //      element on the diagonal.  The columns are not physically interchanged, only relabeled:
  78                 //      indxc[i], the column of the ith pivot element, is the ith column that is reduced, while
  79                 //      indxr[i] is the row in which that pivot element was originally located.  If indxr[i] !=
  80                 //      indxc[i] there is an implied column interchange.  With this form of bookkeeping, the
  81                 //      solution b's will end up in the correct order, and the inverse matrix will be scrambled
  82                 //      by columns.
  83
  84                 if (irow != icol) {
  85                         for (l=1; l<=n; l++)
  86                                 SWAP(a[irow][l], a[icol][l]);
  87                         for (l=1; l<=m; l++)
  88                                 SWAP(b[irow][l], b[icol][l]);
  89                 }
  90
  91                 indxr[i] = irow;
  92                 indxc[i] = icol;
  93                 if (a[icol][icol] == 0) {
  94                         mprintf((0,"Error: Singular matrix-2\n"));
  95                         Int3();
  96                 }
  97                 pivinv = fixdiv(F1_0, a[icol][icol]);
  98                 a[icol][icol] = F1_0;
  99
 100                 for (l=1; l<=n; l++)
 101                         a[icol][l] = fixmul(a[icol][l], pivinv);
 102                 for (l=1; l<=m; l++)
 103                         b[icol][l] = fixmul(b[icol][l], pivinv);
 104
 105                 for (ll=1; ll<=n; ll++)
 106                         if (ll != icol) {
 107                                 dum = a[ll][icol];
 108                                 a[ll][icol] = 0;
 109                                 for (l=1; l<=n; l++)
 110                                         a[ll][l] -= a[icol][l]*dum;
 111                                 for (l=1; l<=m; l++)
 112                                         b[ll][l] -= b[icol][l]*dum;
 113                         }
 114         }
 115
 116         //      This is the end of the main loop over columns of the reduction.  It only remains to unscramble
 117         //      the solution in view of the column interchanges.  We do this by interchanging pairs of
 118         //      columns in the reverse order that the permutation was built up.
 119         for (l=n; l>=1; l--) {
 120                 if (indxr[l] != indxc[l])
 121                         for (k=1; k<=n; k++)
 122                                 SWAP(a[k][indxr[l]], a[k][indxc[l]]);
 123         }
 124
 125 }
 126
 127
 128 //      -----------------------------------------------------------------------------------------------------------------
 129 //      Return true if side is planar, else return false.
 130 int side_is_planar_p(segment *sp, int side)
 131 {
 132         sbyte                   *vp;
 133         vms_vector      *v0,*v1,*v2,*v3;
 134         vms_vector      va,vb;
 135
 136         vp = Side_to_verts[side];
 137         v0 = &Vertices[sp->verts[vp[0]]];
 138         v1 = &Vertices[sp->verts[vp[1]]];
 139         v2 = &Vertices[sp->verts[vp[2]]];
 140         v3 = &Vertices[sp->verts[vp[3]]];
 141
 142         vm_vec_normalize(vm_vec_normal(&va,v0,v1,v2));
 143         vm_vec_normalize(vm_vec_normal(&vb,v0,v2,v3));
 144
 145         // If the two vectors are very close to being the same, then generate one quad, else generate two triangles.
 146         return (vm_vec_dist(&va,&vb) < F1_0/1000);
 147 }
 148
 149 //      -------------------------------------------------------------------------------------------------
 150 //      Return coordinates of a vertex which is vertex v moved so that all sides of which it is a part become planar.
 151 void compute_planar_vert(segment *sp, int side, int v, vms_vector *vp)
 152 {
 153         if ((sp) && (side > -3))
 154                 *vp = Vertices[v];
 155 }
 156
 157 //      -------------------------------------------------------------------------------------------------
 158 //      Making Cursegp:Curside planar.
 159 //      If already planar, return.
 160 //      for each vertex v on side, not part of another segment
 161 //              choose the vertex v which can be moved to make all sides of which it is a part planar, minimizing distance moved
 162 //      if there is no vertex v on side, not part of another segment, give up in disgust
 163 //      Return value:
 164 //              0       curside made planar (or already was)
 165 //              1       did not make curside planar
 166 int make_curside_planar(void)
 167 {
 168         int                     v;
 169         sbyte                   *vp;
 170         vms_vector      planar_verts[4];                        // store coordinates of up to 4 vertices which will make Curside planar, corresponding to each of 4 vertices on side
 171         int                     present_verts[4];                       //      set to 1 if vertex is present
 172
 173         if (side_is_planar_p(Cursegp, Curside))
 174                 return 0;
 175
 176         //      Look at all vertices in side to find a free one.
 177         for (v=0; v<4; v++)
 178                 present_verts[v] = 0;
 179
 180         vp = Side_to_verts[Curside];
 181
 182         for (v=0; v<4; v++) {
 183                 int v1 = vp[v];         // absolute vertex id
 184                 if (is_free_vertex(Cursegp->verts[v1])) {
 185                         compute_planar_vert(Cursegp, Curside, Cursegp->verts[v1], &planar_verts[v]);
 186                         present_verts[v] = 1;
 187                 }
 188         }
 189
 190         //      Now, for each v for which present_verts[v] == 1, there is a vector (point) in planar_verts[v].
 191         //      See which one is closest to the plane defined by the other three points.
 192         //      Nah...just use the first one we find.
 193         for (v=0; v<4; v++)
 194                 if (present_verts[v]) {
 195                         med_set_vertex(vp[v],&planar_verts[v]);
 196                         validate_segment(Cursegp);
 197                         // -- should propagate tmaps to segments or something here...
 198
 199                         return 0;
 200                 }
 201         //      We tried, but we failed, to make Curside planer.
 202         return 1;
 203 }
 204