1
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
2
|
// Copyright 2023 Leszek Koltunski //
|
3
|
// //
|
4
|
// This file is part of Magic Cube. //
|
5
|
// //
|
6
|
// Magic Cube is proprietary software licensed under an EULA which you should have received //
|
7
|
// along with the code. If not, check https://distorted.org/magic/License-Magic-Cube.html //
|
8
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
9
|
|
10
|
package org.distorted.objectlib.helpers;
|
11
|
|
12
|
import org.distorted.library.type.Static3D;
|
13
|
|
14
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
15
|
|
16
|
public class FactoryShape
|
17
|
{
|
18
|
private static final float MAXERROR = 0.0001f;
|
19
|
|
20
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
21
|
// return arithmetic average of the vertices.
|
22
|
|
23
|
private static float[] computeCenterOfMass(float[][] vertices, int num, int[] indices)
|
24
|
{
|
25
|
float[] ret = new float[3];
|
26
|
|
27
|
for( int i=0; i<num; i++ )
|
28
|
{
|
29
|
int ind = indices[i];
|
30
|
float[] v = vertices[ind];
|
31
|
|
32
|
ret[0] += v[0];
|
33
|
ret[1] += v[1];
|
34
|
ret[2] += v[2];
|
35
|
}
|
36
|
|
37
|
ret[0] /= num;
|
38
|
ret[1] /= num;
|
39
|
ret[2] /= num;
|
40
|
|
41
|
return ret;
|
42
|
}
|
43
|
|
44
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
45
|
// Return true iff 'vertex' lies in the plane 'plane'
|
46
|
|
47
|
private static boolean vertexOnPlane(float[] vertex, float[] plane)
|
48
|
{
|
49
|
float d = vertex[0]*plane[0] + vertex[1]*plane[1] + vertex[2]*plane[2] - plane[3];
|
50
|
return (d>=-MAXERROR && d<=MAXERROR);
|
51
|
}
|
52
|
|
53
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
54
|
// return true iff points p0 and p1 are on different sides of plane 'plane'. [ i.e. every way from
|
55
|
// p0 to p1, always has to pass through at least one point from 'plane' ]
|
56
|
|
57
|
private static boolean isOnDifferentSides(float[] p0, float[] p1, float[] plane)
|
58
|
{
|
59
|
float d0 = p0[0]*plane[0] + p0[1]*plane[1] + p0[2]*plane[2] - plane[3];
|
60
|
float d1 = p1[0]*plane[0] + p1[1]*plane[1] + p1[2]*plane[2] - plane[3];
|
61
|
|
62
|
if( d0>=-MAXERROR && d0<=MAXERROR ) return false;
|
63
|
if( d0>MAXERROR ) return d1<-MAXERROR;
|
64
|
return d1>MAXERROR;
|
65
|
}
|
66
|
|
67
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
68
|
// return indices 'indices', but cut to only 'num' of them [original 'indices' array can be longer]
|
69
|
// and arranged so that the 'vertices' those indices are indices for are in counter-clockwise order
|
70
|
// [starting from an arbitrary one].
|
71
|
//
|
72
|
// a) for every vertex, compute the angle the vector from 'center' to this vertex forms with an
|
73
|
// arbitrary line [for example the vector from center to the first vertex]
|
74
|
// b) sort the vertices by this angle [ascending]
|
75
|
|
76
|
private static int[] arrangeCCW(float[] plane, float[] internal, float[][] vertices, int num, int[] indices)
|
77
|
{
|
78
|
float[] angles = new float[num];
|
79
|
float[] center = computeCenterOfMass(vertices,num,indices);
|
80
|
|
81
|
int j= indices[0];
|
82
|
float vx = vertices[j][0] - center[0];
|
83
|
float vy = vertices[j][1] - center[1];
|
84
|
float vz = vertices[j][2] - center[2];
|
85
|
|
86
|
float len = (float)Math.sqrt(vx*vx + vy*vy + vz*vz);
|
87
|
|
88
|
vx /= len;
|
89
|
vy /= len;
|
90
|
vz /= len;
|
91
|
|
92
|
for(int i=0; i<num; i++)
|
93
|
{
|
94
|
int index = indices[i];
|
95
|
float[] v = vertices[index];
|
96
|
|
97
|
float wx = v[0]-center[0];
|
98
|
float wy = v[1]-center[1];
|
99
|
float wz = v[2]-center[2];
|
100
|
|
101
|
float l = (float)Math.sqrt(wx*wx + wy*wy + wz*wz);
|
102
|
|
103
|
wx /= l;
|
104
|
wy /= l;
|
105
|
wz /= l;
|
106
|
|
107
|
angles[i] = computeAngle(plane,internal,vx,vy,vz,wx,wy,wz);
|
108
|
}
|
109
|
|
110
|
int[] ret = new int[num];
|
111
|
for(int i=0; i<num; i++) ret[i] = indices[i];
|
112
|
|
113
|
bubbleSort(angles,ret,num);
|
114
|
|
115
|
return ret;
|
116
|
}
|
117
|
|
118
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
119
|
// strange version of bubbleSort where we sort array 'arr2' according to values of different array
|
120
|
// 'arr1'.
|
121
|
|
122
|
private static void bubbleSort(float[] arr1, int[] arr2, int n)
|
123
|
{
|
124
|
int i, j, temp;
|
125
|
boolean swapped;
|
126
|
float tmp;
|
127
|
|
128
|
for( i=0; i<n-1; i++ )
|
129
|
{
|
130
|
swapped = false;
|
131
|
|
132
|
for( j=0; j<n-i-1; j++ )
|
133
|
{
|
134
|
if( arr1[j] > arr1[j+1] )
|
135
|
{
|
136
|
tmp = arr1[j];
|
137
|
arr1[j] = arr1[j+1];
|
138
|
arr1[j+1]= tmp;
|
139
|
|
140
|
temp = arr2[j];
|
141
|
arr2[j] = arr2[j+1];
|
142
|
arr2[j+1]= temp;
|
143
|
|
144
|
swapped = true;
|
145
|
}
|
146
|
}
|
147
|
|
148
|
if( !swapped ) break;
|
149
|
}
|
150
|
}
|
151
|
|
152
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
153
|
// given a plane with normal vector 'normal', and two in-plane vectors p=(px,py,pz) and
|
154
|
// r=(rx,ry,rz) [in-plane, i.e. cross(p,r) is parallel to 'normal' ] compute the angle between p and r.
|
155
|
//
|
156
|
// the angle belongs to <0,2PI) and goes from p to r counter-clockwise (according to the 'normal' -
|
157
|
// i.e. look at the situation so that the normal points at you and track the angle from p to r)
|
158
|
//
|
159
|
// vectors p and r are guaranteed to be normalized.
|
160
|
//
|
161
|
// 1. compute scalar product of p and r
|
162
|
// 2. result is the cosine of the 'smaller' angle between the two
|
163
|
// 3. a = compute arc cos of this
|
164
|
// 4. return a or (2PI-a), depending on if cross(p,r) agrees in direction with the normal, or not.
|
165
|
// 'point' is the 'internalPoint' --> the normal needs to point away from this, which is always
|
166
|
// the case in case of the external walls, but not in case of the internal cuts.
|
167
|
|
168
|
private static float computeAngle(float[] plane, float[] point,
|
169
|
float px, float py, float pz, float rx, float ry, float rz)
|
170
|
{
|
171
|
float s = px*rx + py*ry + pz*rz;
|
172
|
|
173
|
if( s> 1 ) s= 1;
|
174
|
if( s<-1 ) s=-1;
|
175
|
|
176
|
float a = (float) Math.acos(s);
|
177
|
float[] cross = crossProduct(px,py,pz,rx,ry,rz);
|
178
|
|
179
|
float scalar1 = cross[0]*plane[0] + cross[1]*plane[1] + cross[2]*plane[2];
|
180
|
float scalar2 = point[0]*plane[0] + point[1]*plane[1] + point[2]*plane[2];
|
181
|
|
182
|
return ((scalar1>=0)^(scalar2>=plane[3])) ? a : (float)(2*Math.PI - a);
|
183
|
}
|
184
|
|
185
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
186
|
|
187
|
private static float[] crossProduct(float a1, float a2, float a3, float b1, float b2, float b3)
|
188
|
{
|
189
|
float[] ret = new float[3];
|
190
|
|
191
|
ret[0] = a2*b3 - a3*b2;
|
192
|
ret[1] = a3*b1 - a1*b3;
|
193
|
ret[2] = a1*b2 - a2*b1;
|
194
|
|
195
|
return ret;
|
196
|
}
|
197
|
|
198
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
199
|
// if normals to planes p1,p2,p3 are linearly independent, than those planes intersect in one point.
|
200
|
// Return this point. Otherwise, return null.
|
201
|
|
202
|
private static float[] computeCommonPoint(float[] p1, float[] p2, float[] p3)
|
203
|
{
|
204
|
float det = determinant(p1[0],p1[1],p1[2],p2[0],p2[1],p2[2],p3[0],p3[1],p3[2]);
|
205
|
|
206
|
if( det>=-MAXERROR && det<=MAXERROR ) return null;
|
207
|
|
208
|
float deX = determinant(p1[3],p1[1],p1[2],p2[3],p2[1],p2[2],p3[3],p3[1],p3[2]);
|
209
|
float deY = determinant(p1[0],p1[3],p1[2],p2[0],p2[3],p2[2],p3[0],p3[3],p3[2]);
|
210
|
float deZ = determinant(p1[0],p1[1],p1[3],p2[0],p2[1],p2[3],p3[0],p3[1],p3[3]);
|
211
|
|
212
|
return new float[] { deX/det, deY/det, deZ/det };
|
213
|
}
|
214
|
|
215
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
216
|
// a1 a2 a3 a1 a2
|
217
|
// b1 b2 b3 b1 b2
|
218
|
// c1 c2 c3 c1 c2
|
219
|
|
220
|
private static float determinant(float a1, float a2, float a3,
|
221
|
float b1, float b2, float b3,
|
222
|
float c1, float c2, float c3 )
|
223
|
{
|
224
|
float det = 0;
|
225
|
|
226
|
det += a1*b2*c3;
|
227
|
det += a2*b3*c1;
|
228
|
det += a3*b1*c2;
|
229
|
det -= a3*b2*c1;
|
230
|
det -= a1*b3*c2;
|
231
|
det -= a2*b1*c3;
|
232
|
|
233
|
return det;
|
234
|
}
|
235
|
|
236
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
237
|
// how many faces does the point (x,y,z) belong to? (i.e. 3-->corner, 2-->edge, 1--> face point,
|
238
|
// 0 --> internal
|
239
|
|
240
|
private static int numberInFaces(Static3D[] faceAxis, float[] dist3D, float x, float y, float z)
|
241
|
{
|
242
|
int numA = faceAxis.length;
|
243
|
int ret = 0;
|
244
|
|
245
|
for(int a=0; a<numA; a++)
|
246
|
{
|
247
|
Static3D ax = faceAxis[a];
|
248
|
float dist = x*ax.get0() + y*ax.get1() + z*ax.get2();
|
249
|
float diff = dist-dist3D[a];
|
250
|
if( diff>-MAXERROR && diff<MAXERROR ) ret++;
|
251
|
}
|
252
|
|
253
|
return ret;
|
254
|
}
|
255
|
|
256
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
257
|
// does the (x,y,z) point belong to any of the faces defined by (faceAxis,dist) ?
|
258
|
|
259
|
private static boolean isInFaces(Static3D[] faceAxis, float[] dist3D, float x, float y, float z)
|
260
|
{
|
261
|
int numA = faceAxis.length;
|
262
|
|
263
|
for(int a=0; a<numA; a++)
|
264
|
{
|
265
|
Static3D ax = faceAxis[a];
|
266
|
float dist = x*ax.get0() + y*ax.get1() + z*ax.get2();
|
267
|
float diff = dist-dist3D[a];
|
268
|
if( diff>-MAXERROR && diff<MAXERROR ) return true;
|
269
|
}
|
270
|
|
271
|
return false;
|
272
|
}
|
273
|
|
274
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
275
|
|
276
|
private static boolean isDifferent(float[][] table, int num, float[] vert)
|
277
|
{
|
278
|
for( int i=0; i<num; i++ )
|
279
|
{
|
280
|
float[] t = table[i];
|
281
|
float dx = t[0]-vert[0];
|
282
|
float dy = t[1]-vert[1];
|
283
|
float dz = t[2]-vert[2];
|
284
|
float dist = dx*dx + dy*dy + dz*dz;
|
285
|
|
286
|
if( dist>=-MAXERROR && dist<=MAXERROR ) return false;
|
287
|
}
|
288
|
|
289
|
return true;
|
290
|
}
|
291
|
|
292
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
293
|
// PUBLIC API
|
294
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
295
|
|
296
|
public static float[][] computePotentialVertices(float[][] cutPlanes)
|
297
|
{
|
298
|
int numPlanes = cutPlanes.length;
|
299
|
int numRet = numPlanes*(numPlanes-1)*(numPlanes-2)/6;
|
300
|
float[][] ret = new float[numRet][];
|
301
|
int index = 0;
|
302
|
|
303
|
for(int i=0; i<numPlanes; i++)
|
304
|
for(int j=i+1; j<numPlanes; j++)
|
305
|
for(int k=j+1; k<numPlanes; k++)
|
306
|
{
|
307
|
float[] vert = computeCommonPoint(cutPlanes[i],cutPlanes[j],cutPlanes[k]);
|
308
|
|
309
|
if( vert!=null && isDifferent(ret,index,vert) )
|
310
|
{
|
311
|
ret[index] = vert;
|
312
|
index++;
|
313
|
}
|
314
|
}
|
315
|
|
316
|
if( index<numRet )
|
317
|
{
|
318
|
float[][] ret2 = new float[index][];
|
319
|
for(int i=0; i<index; i++) ret2[i] = ret[i];
|
320
|
return ret2;
|
321
|
}
|
322
|
|
323
|
return ret;
|
324
|
}
|
325
|
|
326
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
327
|
// return a table of indices.length ints.
|
328
|
// Nth int is 0 if the face defined by indices[N] is an external one; 1 otherwise.
|
329
|
|
330
|
public static int[] produceBandIndices(float[][] vertices, float[] position, int[][] indices, Static3D[] faceAxis, float[] dist3D )
|
331
|
{
|
332
|
int numF = indices.length;
|
333
|
int[] ret = new int[numF];
|
334
|
float x,y,z;
|
335
|
|
336
|
for(int f=0; f<numF; f++)
|
337
|
{
|
338
|
int[] faceI = indices[f];
|
339
|
x = y = z = 0;
|
340
|
|
341
|
for(int i : faceI)
|
342
|
{
|
343
|
float[] vertex = vertices[i];
|
344
|
x += vertex[0];
|
345
|
y += vertex[1];
|
346
|
z += vertex[2];
|
347
|
}
|
348
|
|
349
|
int numV = faceI.length;
|
350
|
x /= numV;
|
351
|
y /= numV;
|
352
|
z /= numV;
|
353
|
|
354
|
x += position[0];
|
355
|
y += position[1];
|
356
|
z += position[2];
|
357
|
|
358
|
ret[f] = isInFaces(faceAxis,dist3D, x,y,z) ? 0 : 1;
|
359
|
}
|
360
|
|
361
|
return ret;
|
362
|
}
|
363
|
|
364
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
365
|
// return a table of vertices.length ints.
|
366
|
// Nth int is 0 if vertices[N] lies in a corner or edge of the shape; -1 otherwise.
|
367
|
|
368
|
public static int[] computeVertexEffectsIndices(float[][] vertices, float[] position, Static3D[] faceAxis, float[] dist3D)
|
369
|
{
|
370
|
int numV = vertices.length;
|
371
|
int[] ret = new int[numV];
|
372
|
|
373
|
for(int v=0; v<numV; v++)
|
374
|
{
|
375
|
float[] vert = vertices[v];
|
376
|
|
377
|
float x = vert[0] + position[0];
|
378
|
float y = vert[1] + position[1];
|
379
|
float z = vert[2] + position[2];
|
380
|
|
381
|
ret[v] = numberInFaces(faceAxis,dist3D,x,y,z)>=2 ? 0 : -1;
|
382
|
}
|
383
|
|
384
|
return ret;
|
385
|
}
|
386
|
|
387
|
///////////////////////////////////////////////////////////////////////////////////////////////////
|
388
|
// 1. for every potentialVertex, if it is on the same side of every plane present as the
|
389
|
// 'internalPoint', it must be one of the vertices of the ObjectShape. Remember which planes it
|
390
|
// belongs to and add it to the list of vertices.
|
391
|
// 2. Having complete list of vertices and a bitmap of which planes each vertex belongs to, for each
|
392
|
// plane:
|
393
|
// a) assemble the list of vertices that belong to this plane
|
394
|
// b) compute their center of mass and using this as a pivot, arrange the list into CCW order
|
395
|
// c) out of the vertices and the CCW'd indices, construct an ObjectShape.
|
396
|
|
397
|
public static ObjectShape createShape(float[][] cutPlanes, float[][] potentialVertices, float[] position, float[] internalPoint)
|
398
|
{
|
399
|
int numV = potentialVertices.length;
|
400
|
float[][] tmp = new float[numV][];
|
401
|
int numVertices = 0;
|
402
|
|
403
|
for( float[] vert : potentialVertices )
|
404
|
{
|
405
|
boolean ok = true;
|
406
|
|
407
|
for(float[] plane : cutPlanes)
|
408
|
{
|
409
|
if( isOnDifferentSides(internalPoint, vert, plane) )
|
410
|
{
|
411
|
ok = false;
|
412
|
break;
|
413
|
}
|
414
|
}
|
415
|
|
416
|
if( ok )
|
417
|
{
|
418
|
tmp[numVertices] = vert;
|
419
|
numVertices++;
|
420
|
}
|
421
|
}
|
422
|
|
423
|
float[][] vertices = new float[numVertices][];
|
424
|
|
425
|
for(int v=0; v<numVertices; v++)
|
426
|
{
|
427
|
float[] t = tmp[v];
|
428
|
vertices[v] = new float[] { t[0],t[1],t[2] };
|
429
|
}
|
430
|
|
431
|
int numP = cutPlanes.length;
|
432
|
int[] numVerticesOnPlane = new int[numP];
|
433
|
int[][] verticesOnPlane = new int[numP][numVertices];
|
434
|
|
435
|
for(int p=0; p<numP; p++)
|
436
|
{
|
437
|
float[] plane = cutPlanes[p];
|
438
|
|
439
|
for(int v=0; v<numVertices; v++)
|
440
|
{
|
441
|
float[] vert = vertices[v];
|
442
|
|
443
|
if( vertexOnPlane(vert,plane) )
|
444
|
{
|
445
|
verticesOnPlane[p][numVerticesOnPlane[p]] = v;
|
446
|
numVerticesOnPlane[p]++;
|
447
|
}
|
448
|
}
|
449
|
}
|
450
|
|
451
|
int numFaces = 0;
|
452
|
for(int p=0; p<numP; p++)
|
453
|
if( numVerticesOnPlane[p]>=3 ) numFaces++;
|
454
|
|
455
|
int[][] indices = new int[numFaces][];
|
456
|
int index = 0;
|
457
|
|
458
|
for(int p=0; p<numP; p++ )
|
459
|
if( numVerticesOnPlane[p]>=3 )
|
460
|
{
|
461
|
indices[index] = arrangeCCW(cutPlanes[p],internalPoint,vertices,numVerticesOnPlane[p],verticesOnPlane[p]);
|
462
|
index++;
|
463
|
}
|
464
|
|
465
|
for( float[] v: vertices )
|
466
|
{
|
467
|
v[0] -= position[0];
|
468
|
v[1] -= position[1];
|
469
|
v[2] -= position[2];
|
470
|
}
|
471
|
|
472
|
return new ObjectShape(vertices,indices);
|
473
|
}
|
474
|
}
|