add code to read MRI data
[anneal.git] / src / edu / berkeley / qfat / Main.java
1 package edu.berkeley.qfat;
2 import java.awt.*;
3 import java.io.*;
4 import java.awt.event.*;
5 import javax.swing.*;
6 import javax.media.opengl.*;
7 import javax.media.opengl.glu.*;
8 import java.util.*;
9 import edu.berkeley.qfat.bind.*;
10 import edu.berkeley.qfat.geom.*;
11 import edu.berkeley.qfat.stl.*;
12 import edu.berkeley.qfat.voxel.*;
13 import edu.berkeley.qfat.geom.Point;
14 import edu.berkeley.qfat.geom.Polygon;
15
16 /*
17
18
19 Todo
20 - review catmull-clark; move vertex points?
21 - re-anneal fish
22
23 - show constraints (?)
24 - show in red if user tries to move a point to an illegal location
25
26 - eliminate use of bindinggroupschanged()
27
28 - post qfat on my software page
29
30
31 With Sequin
32 - constraints admit intersections (lattice)
33 - constraints may be transformed linearly
34
35
36
37
38 Log console
39   - blend-shaded overlay?  slick.
40
41 face/vertex count
42 rendering FPS
43 ability to not draw edges between faces
44
45
46 three circumcircles showing crystal ball -- these don't get scaled
47 axes?
48
49 drawing modes:
50   - bounding box
51   - vertices
52   - edges
53   - visible-edges
54   - flat with or without edges
55   - shaded with or without edges
56   * contrasting-faces
57
58
59 quadric decimation?
60
61 show normals
62 show bounding box
63 show axes (big+fat)
64  */
65
66 // TO DO:
67 //
68 // - Ability to snap three views to orthgonal
69 // - SLIDE UI
70 //     - left button -> crystal ball
71 //     - translate
72 //     - rightbutton/mousewheel zoom
73 // - v+click to select vertex
74 // - show: constraints, bindings, faces
75 //
76 // Editing:
77 //  - fracture edge, face
78 //  - change aspect ratio
79 //  - translate, rotate goal mesh
80 //  - ability to select a point, rotate the model, then move the point
81 //  - when moving a vertex in one window, show that window's axes in all other windows
82 //
83
84 /*
85   blender keys
86   - middle mouse = option+click
87   - right mouse = command+click
88
89   3,7,1 = view along axes (control for opp direction)
90   4, 8, 7, 2 = rotate in discrete increments (+control to translate)
91   middle trag: rotate space
92   shift+middle drag: translate space
93   wheel: zoom
94   home: home view: take current angle, zoom to whole scnee
95   5 = ortho vs non-ortho
96 */
97
98 /*
99 Meshlab Notes:
100 Log console
101   - blend-shaded overlay?  slick.
102
103 face/vertex count
104 rendering FPS
105 ability to not draw edges between faces
106
107
108 three circumcircles showing crystal ball -- these don't get scaled
109 axes?
110
111 drawing modes:
112   - bounding box
113   - vertices
114   - edges
115   - visible-edges
116   - flat with or without edges
117   - shaded with or without edges
118   * contrasting-faces
119
120
121 quadric decimation?
122
123 show normals
124 show bounding box
125 show axes (big+fat)
126  */
127
128 public class Main extends InteractiveMeshViewer {
129
130     public static int verts = 1;
131
132     public static final Random random = new Random();
133     
134     /** magnification factor */
135     private static final float MAG = 1;
136     public static final float MATCHING_EPSILON = 0.001f;
137
138     private static boolean small(float f) { return Math.abs(f) < 0.001; }
139     public void generateTile(Matrix[] matrices, Mesh mesh) {
140         mesh.coalesce = true;
141         HashSet<HalfSpace> halfSpaces = new HashSet<HalfSpace>();
142         HashSet<Polygon> polygons = new HashSet<Polygon>();
143         for(Matrix m : matrices) {
144                 Vec v = m.getTranslationalComponent();
145                 if (v.mag() < 0.0001) continue;
146                 v = v.times(-1);
147                 v = v.times(0.5f);
148                 Point p = Point.ZERO.plus(v);
149                 v = v.times(-1);
150                 
151                 //System.out.println(v);
152                 HalfSpace hs = new HalfSpace(p, v.norm());
153                 halfSpaces.add(hs);
154                 polygons.add(new Polygon(hs));
155         }
156         for(Polygon p : polygons) {
157             System.out.println(p.plane.norm() + " " + p.plane.d);
158             for(HalfSpace hs : halfSpaces) {
159                 if (p.plane==hs) continue;
160                 p = p.intersect(hs);
161             }
162             p.tesselate(mesh);
163         }
164     }
165
166     private void quad(Mesh mesh, Matrix m, Point p1_, Point p2_, Point p3_, Point p4_) {
167         Point p1 = m.times(p1_);
168         Point p2 = m.times(p2_);
169         Point p3 = m.times(p3_);
170         Point p4 = m.times(p4_);
171         Point c  = new Point((p1.x+p2.x+p3.x+p4.x)/4,
172                              (p1.y+p2.y+p3.y+p4.y)/4,
173                              (p1.z+p2.z+p3.z+p4.z)/4);
174         mesh.newT(p1, p2, c, null, 0);
175         mesh.newT(p2, p3, c, null, 0);
176         mesh.newT(p3, p4, c, null, 0);
177         mesh.newT(p4, p1, c, null, 0);
178     }
179
180     public void loadGoal(String file) {
181         try {
182             StlFile stlf = new StlFile();
183             InputStream res = this.getClass().getClassLoader().getResourceAsStream(file);
184             stlf.readBinaryFile(file, res);
185             setGoal(new Mesh(false));
186             for(int i=0; i<stlf.coordArray.length; i+=3) {
187                 Point p0 = new Point(stlf.coordArray[i+0].x * MAG, stlf.coordArray[i+0].y * MAG, stlf.coordArray[i+0].z * MAG);
188                 Point p1 = new Point(stlf.coordArray[i+1].x * MAG, stlf.coordArray[i+1].y * MAG, stlf.coordArray[i+1].z * MAG);
189                 Point p2 = new Point(stlf.coordArray[i+2].x * MAG, stlf.coordArray[i+2].y * MAG, stlf.coordArray[i+2].z * MAG);
190                 Vec n    = new Vec(stlf.normArray[i/3].x * MAG, stlf.normArray[i/3].y  * MAG, stlf.normArray[i/3].z * MAG);
191                 Mesh.T t  = goal.newT(p0, p1, p2, n, 0);
192             }
193             float goal_width  = goal.diagonal().dot(new Vec(1, 0, 0));
194             float goal_height = goal.diagonal().dot(new Vec(0, 1, 0));
195             float goal_depth  = goal.diagonal().dot(new Vec(0, 0, 1));
196             
197         } catch (Exception e) { throw new RuntimeException(e);}
198     }
199     public void fixupGoal() { fixupGoal(true, true); }
200     public void fixupGoal(boolean recenter, boolean lock) {
201         // translate to match centroid
202         if (recenter)
203             goal.transform(Matrix.translate(tile.centroid().minus(goal.centroid())));
204         if (lock) {
205             goal.makeVerticesImmutable();
206             tile.error_against = goal;
207             goal.error_against = tile;
208         }
209     }
210
211     public Main(JFrame f) { super(f); }
212
213     public synchronized void fixupTile() { 
214         for(Matrix m1 : transforms) {
215             for(Matrix m2 : transforms) {
216                 if (m1==m2) continue;
217                 for(Mesh.T t1 : tile) {
218                     for(Mesh.T t2 : tile) {
219
220                         Matrix m = m1.inverse().times(m2);
221                         Point t1v1 = m.times(t1.v1().p);
222                         Point t1v2 = m.times(t1.v2().p);
223                         Point t1v3 = m.times(t1.v3().p);
224
225                         if (t1v1.distance(t2.v1().p) < MATCHING_EPSILON &&
226                             t1v2.distance(t2.v3().p) < MATCHING_EPSILON &&
227                             t1v3.distance(t2.v2().p) < MATCHING_EPSILON) {
228                             t2.e3().bindEdge(t1.e1().pair, m);
229                             t2.e2().bindEdge(t1.e2().pair, m);
230                             t2.e1().bindEdge(t1.e3().pair, m);
231                         }
232                         if (t1v2.distance(t2.v1().p) < MATCHING_EPSILON &&
233                             t1v3.distance(t2.v3().p) < MATCHING_EPSILON &&
234                             t1v1.distance(t2.v2().p) < MATCHING_EPSILON) {
235                             t2.e3().bindEdge(t1.e2().pair, m);
236                             t2.e2().bindEdge(t1.e3().pair, m);
237                             t2.e1().bindEdge(t1.e1().pair, m);
238                         }
239                         if (t1v3.distance(t2.v1().p) < MATCHING_EPSILON &&
240                             t1v1.distance(t2.v3().p) < MATCHING_EPSILON &&
241                             t1v2.distance(t2.v2().p) < MATCHING_EPSILON) {
242                             t2.e3().bindEdge(t1.e3().pair, m);
243                             t2.e2().bindEdge(t1.e1().pair, m);
244                             t2.e1().bindEdge(t1.e2().pair, m);
245                         }
246
247                         if (t1v1.distance(t2.v1().p) < MATCHING_EPSILON &&
248                             t1v2.distance(t2.v2().p) < MATCHING_EPSILON &&
249                             t1v3.distance(t2.v3().p) < MATCHING_EPSILON) {
250                             t2.e1().bindEdge(t1.e1(), m);
251                             t2.e2().bindEdge(t1.e2(), m);
252                             t2.e3().bindEdge(t1.e3(), m);
253                         }
254                         if (t1v2.distance(t2.v1().p) < MATCHING_EPSILON &&
255                             t1v3.distance(t2.v2().p) < MATCHING_EPSILON &&
256                             t1v1.distance(t2.v3().p) < MATCHING_EPSILON) {
257                             t2.e2().bindEdge(t1.e1(), m);
258                             t2.e3().bindEdge(t1.e2(), m);
259                             t2.e1().bindEdge(t1.e3(), m);
260                         }
261                         if (t1v3.distance(t2.v1().p) < MATCHING_EPSILON &&
262                             t1v1.distance(t2.v2().p) < MATCHING_EPSILON &&
263                             t1v2.distance(t2.v3().p) < MATCHING_EPSILON) {
264                             t2.e3().bindEdge(t1.e1(), m);
265                             t2.e1().bindEdge(t1.e2(), m);
266                             t2.e2().bindEdge(t1.e3(), m);
267                         }
268
269                     }
270                 }
271             }
272         }
273         tile.rebindPoints();
274         System.out.println("tile volume: " + tile.volume());
275         System.out.println("goal volume: " + goal.volume());
276         tile.error_against = goal;
277         goal.error_against = tile;
278         fixupGoal();
279     }
280
281
282     public class MyMenuItem extends JMenuItem implements ActionListener {
283         public MyMenuItem(String s) {
284             super(s);
285             this.addActionListener(this);
286         }
287         public void actionPerformed(ActionEvent event) {
288             synchronized(Main.this) {
289                     hit();
290             }
291         }
292         public void hit() {}
293     }
294
295     public void marchingCubes() {
296         Mesh mesh = new Mesh(false);
297         mesh.coalesce = true;
298         MarchingCubes.march(new VoxelData() {
299                 float radius = 1.0f;
300                 public float getMaxX()        { return  1.0f; }
301                 public float getMinX()        { return -1.0f; }
302                 public int   getNumSamplesX() { return   10; }
303                 public float getMaxY()        { return  1.0f; }
304                 public float getMinY()        { return -1.0f; }
305                 public int   getNumSamplesY() { return   10; }
306                 public float getMaxZ()        { return  1.0f; }
307                 public float getMinZ()        { return -1.0f; }
308                 public int   getNumSamplesZ() { return   10; }
309                 public float getSample(Point p) {
310                     double x = p.x;
311                     double y = p.y;
312                     double z = p.z;
313                     return (float)(radius-Math.sqrt(x*x+y*y+z*z));
314                 }
315             },
316             mesh);
317         setTile(mesh);
318         //fixupTile();
319     }
320
321     public void marchingCubes2() {
322         try {
323
324             final float[][][] samples = new float[256][256][124];
325             double total = 0;
326             int count = 0;
327             for(int i=1; i<=124; i++) {
328                 String ix = ""+i;
329                 while(ix.length() < 3) ix = "0"+ix;
330                 FileInputStream fis = new FileInputStream("/Users/megacz/Desktop/projects/mri brain images/spgr/I."+ix);
331                 DataInputStream dis = new DataInputStream(new BufferedInputStream(fis));
332                 dis.skip(7904);
333                 for(int x=0; x<256; x++)
334                     for(int y=0; y<256; y++) {
335                         short s = dis.readShort();
336                         samples[x][y][i-1] = (float)s;
337                         total += samples[x][y][i-1];
338                         count++;
339                     }
340             }
341             System.out.println("done reading samples; average sample is " + (total/count));
342
343             /*
344             PrintWriter pw = new PrintWriter(new FileOutputStream("/tmp/out"));
345             pw.println("new short[][][] {");
346             for(int x=0; x<256; x++) {
347                 if (x>0) pw.println("  ,");
348                 pw.println("  {");
349                 for(int y=0; y<256; y++) {
350                     if (y>0) pw.println("    ,");
351                     pw.print("    {");
352                     for(int z=0; z<124; z++) {
353                         if (z>0) pw.print(", ");
354                         if ((z % 20) == 0) {
355                             pw.println();
356                             pw.print("      ");
357                         }
358                         pw.print(shorts[x][y][z]);
359                     }
360                     pw.println();
361                     pw.println("    }");
362                 }
363                 pw.println("  }");
364             }
365             pw.println("};");
366             pw.flush();
367             pw.close();
368             */
369
370             Mesh mesh = new Mesh(false);
371             mesh.coalesce = true;
372             MarchingCubes.march(new VoxelData() {
373                     float radius = 1.0f;
374                     public float getMaxX()        { return  2.4f; }
375                     public float getMinX()        { return    0f; }
376                     public int   getNumSamplesX() { return   256/10; }
377                     public float getMaxY()        { return  2.4f; }
378                     public float getMinY()        { return    0f; }
379                     public int   getNumSamplesY() { return   256/10; }
380                     public float getMaxZ()        { return 1.86f; }
381                     public float getMinZ()        { return    0f; }
382                     public int   getNumSamplesZ() { return   124/10; }
383                     public float getSample(Point p) {
384                         int x = (int)Math.floor(p.x / ((double)getMaxX() / (double)256));
385                         int y = (int)Math.floor(p.y / ((double)getMaxY() / (double)256));
386                         int z = (int)Math.floor(p.z / ((double)getMaxZ() / (double)124));
387                         if ( (x<0) || (x>=samples.length) ) return 0;
388                         if ( (y<0) || (y>=samples[x].length) ) return 0;
389                         if ( (z<0) || (z>=samples[x][y].length) ) return 0;
390                         return samples[x][y][z];
391                     }
392                 },
393                 200,
394                 mesh);
395             setTile(mesh);
396             //fixupTile();
397
398         } catch (Exception e) {
399             e.printStackTrace();
400         }
401     }
402
403     public void hexBrick(boolean offset, boolean rotated) {
404         setTile(new Mesh(false));
405                 float width  = (float)0.8;
406                 float depth  = (float)0.08;
407                 float height = (float)0.4;
408                 float rshift =   width/2;
409                 float lshift = -(width/2);
410                 float halfup = 0;
411                 Point ltf = new Point(lshift,  (depth/2),  (height/2));
412                 Point mtf = new Point( 0.0,    (depth/2),  (height/2));
413                 Point rtf = new Point(rshift,  (depth/2),  (height/2));
414                 Point lbf = new Point(lshift, -(depth/2),  (height/2));
415                 Point mbf = new Point( 0.0,   -(depth/2),  (height/2));
416                 Point rbf = new Point(rshift, -(depth/2),  (height/2));
417
418                 Point ltc = new Point(lshift,  (depth/2), 0);
419                 Point mtc = new Point( 0.0,    (depth/2), 0);
420                 Point rtc = new Point(rshift,  (depth/2), 0);
421                 Point lbc = new Point(lshift, -(depth/2), 0);
422                 Point mbc = new Point( 0.0,   -(depth/2), 0);
423                 Point rbc = new Point(rshift, -(depth/2), 0);
424
425                 Point ltn = new Point(lshift,  (depth/2), -(height/2));
426                 Point mtn = new Point( 0.0,    (depth/2), -(height/2));
427                 Point rtn = new Point(rshift,  (depth/2), -(height/2));
428                 Point lbn = new Point(lshift, -(depth/2), -(height/2));
429                 Point mbn = new Point( 0.0,   -(depth/2), -(height/2));
430                 Point rbn = new Point(rshift, -(depth/2), -(height/2));
431
432         
433                 Point[] points = new Point[] {
434                     ltf,
435                     mtf,
436                     rtf,
437                     lbf,
438                     mbf,
439                     rbf,
440
441                     ltc,
442                     mtc,
443                     rtc,
444                     lbc,
445                     mbc,
446                     rbc,
447
448                     ltn,
449                     mtn,
450                     rtn,
451                     lbn,
452                     mbn,
453                     rbn
454                 };
455
456
457                 // top
458                 tile.newT(ltf, mtf, mtc, null, 1);
459                 tile.newT(mtc, ltc, ltf, null, 1);
460                 tile.newT(mtf, rtf, rtc, null, 1);
461                 tile.newT(rtc, mtc, mtf, null, 1);
462
463                 tile.newT(ltc, mtc, mtn, null, 1);
464                 tile.newT(mtn, ltn, ltc, null, 1);
465                 tile.newT(mtc, rtc, rtn, null, 1);
466                 tile.newT(rtn, mtn, mtc, null, 1);
467
468                 // bottom (swap normals)
469                 tile.newT(mbf, lbf, mbc, null, 2);
470                 tile.newT(lbc, mbc, lbf, null, 2);
471                 tile.newT(rbf, mbf, rbc, null, 2);
472                 tile.newT(mbc, rbc, mbf, null, 2);
473
474                 tile.newT(mbc, lbc, mbn, null, 2);
475                 tile.newT(lbn, mbn, lbc, null, 2);
476
477                 tile.newT(rbc, mbc, rbn, null, 2);
478                 tile.newT(mbn, rbn, mbc, null, 2);
479
480
481                 // left
482                 tile.newT(ltf, ltc, lbc, null, 3);
483                 tile.newT(lbc, lbf, ltf, null, 3);
484                 tile.newT(ltc, ltn, lbn, null, 3);
485                 tile.newT(lbn, lbc, ltc, null, 3);
486
487                 // right (swap normals)
488                 tile.newT(rtc, rtf, rbc, null, 4);
489                 tile.newT(rbf, rbc, rtf, null, 4);
490                 tile.newT(rtn, rtc, rbn, null, 4);
491                 tile.newT(rbc, rbn, rtc, null, 4);
492
493                 // front
494                 tile.newT(ltn, mtn, mbn, null, 5);
495                 tile.newT(ltn, mbn, lbn, null, 5);
496                 tile.newT(mtn, rtn, rbn, null, 5);
497                 tile.newT(mtn, rbn, mbn, null, 5);
498
499                 // back
500                 tile.newT(mtf, ltf, mbf, null, 6);
501                 tile.newT(mbf, ltf, lbf, null, 6);
502                 tile.newT(rtf, mtf, rbf, null, 6);
503                 tile.newT(rbf, mtf, mbf, null, 6);
504
505                 if (offset) {
506                   transforms = new Matrix[] {
507                       Matrix.translate(new Vec(lshift,      0,    height)),
508                       Matrix.translate(new Vec(lshift,      0,   -height)),
509                       Matrix.translate(new Vec(rshift,      0,    height)),
510                       Matrix.translate(new Vec(rshift,      0,   -height)),
511                       Matrix.translate(new Vec( width,      0,         0)),
512                       Matrix.translate(new Vec(-width,      0,         0)),
513                       Matrix.translate(new Vec(lshift,  depth, 0)),
514                       Matrix.translate(new Vec(lshift, -depth, 0)),
515                       Matrix.translate(new Vec(rshift,  depth, 0)),
516                       Matrix.translate(new Vec(rshift, -depth, 0)),
517                       Matrix.ONE,
518                   };
519                 } else if (rotated) {
520                     HashSet<Mesh.E> es = new HashSet<Mesh.E>();
521                     for(Mesh.T t : tile) {
522                         es.add(t.e1());
523                         es.add(t.e2());
524                         es.add(t.e3());
525                     }
526                     for(Mesh.E e : es) {
527                         if (e.v1.p.x == e.v2.p.x && e.v1.p.y == e.v2.p.y) continue;
528                         if (e.v1.p.z == e.v2.p.z && e.v1.p.y == e.v2.p.y) continue;
529                         if (e.v1.p.x == e.v2.p.x && e.v1.p.z == e.v2.p.z) continue;
530                         e.shatter();
531                     }
532                     transforms = new Matrix[] {
533                         Matrix.translate(new Vec(0,   0,    height)).times(Matrix.rotate(new Vec(0, 0, 1), (float)Math.PI)),
534                         Matrix.translate(new Vec(0,   0,   -height)).times(Matrix.rotate(new Vec(0, 0, 1), (float)Math.PI)),
535                         Matrix.translate(new Vec(0,  depth, 0)),
536                         Matrix.translate(new Vec(0, -depth, 0)),
537                         Matrix.translate(new Vec( width,  0,    0)),
538                         Matrix.translate(new Vec(-width,  0,    0)),
539                       Matrix.ONE,
540                   };
541                 } else {
542                   transforms = new Matrix[] {
543                       Matrix.translate(new Vec(lshift,      0,    height)),
544                       Matrix.translate(new Vec(lshift,      0,   -height)),
545                       Matrix.translate(new Vec(rshift,      0,    height)),
546                       Matrix.translate(new Vec(rshift,      0,   -height)),
547                       Matrix.translate(new Vec( width,      0,         0)),
548                       Matrix.translate(new Vec(-width,      0,         0)),
549                       Matrix.translate(new Vec(0,  depth, 0)),
550                       Matrix.translate(new Vec(0, -depth, 0)),
551                       Matrix.ONE,
552                   };
553                 }
554                 
555                 fixupTile();
556     }
557
558     public class MyMenuBar extends JMenuBar {
559
560         public MyMenuBar() {
561
562             JMenu tileMenu = new JMenu("Tile");
563             JMenu goalMenu = new JMenu("Goal");
564             JMenu hideMenu = new JMenu("Actions");
565
566             hideMenu.add(new MyMenuItem("Start Anneal") { public void hit() { anneal = true; }});
567             hideMenu.add(new MyMenuItem("Stop Anneal") { public void hit() { anneal = false; }});
568             hideMenu.add(new MyMenuItem("Reset to high temperature") { public void hit() { temp = 1; }});
569             hideMenu.add(new MyMenuItem("Subdivide surface") { public void hit() { breaks++; }});
570             hideMenu.add(new MyMenuItem("Show Goal") { public void hit() { goalon = true; }});
571             hideMenu.add(new MyMenuItem("Hide Goal") { public void hit() { goalon = false; }});
572             hideMenu.add(new MyMenuItem("Show All Neighbors") { public void hit() { neighbors = true; }});
573             hideMenu.add(new MyMenuItem("Show One Neighbor Wireframe") { public void hit() { neighborsWireOne = true; }});
574             hideMenu.add(new MyMenuItem("Show All Neighbors Wireframe") { public void hit() { neighborsWire = true; neighborsWireOne = false;}});
575             hideMenu.add(new MyMenuItem("Hide Neighbors") { public void hit() { neighborsWire = false; neighborsWireOne = false; neighbors = false; }});
576
577             goalMenu.add(new MyMenuItem("Fish with face") { public void hit() {
578                 loadGoal("face.stl");
579                 //goal.transform(Matrix.rotate(new Vec(0, 1, 0), (float)(Math.PI/2)));
580                 goal.transform(Matrix.rotate(new Vec(0, 0, 1), (float)(Math.PI/2)));
581                 //goal.transform(Matrix.scale(1, 2.2f, 1));
582                 float factor = (float)Math.pow(tile.volume() / goal.volume(), 1.0/3.0);
583                 //factor = factor * 0.8f;
584                 goal.transform(Matrix.scale(0.3f));
585                 //goal.transform(Matrix.rotate(new Vec(0,1,0), (float)(Math.PI*1/3)));
586                 goal.transform(Matrix.rotate(new Vec(1,0,0), (float)(Math.PI*1/3)));
587                 fixupGoal(true, false);
588                 //goal.transform(Matrix.translate(new Vec(0.145f, 0, 0)));
589                 fixupGoal(false, true);
590             }});
591             goalMenu.add(new MyMenuItem("Fish") { public void hit() {
592                 loadGoal("fish.stl");
593                 //goal.transform(Matrix.rotate(new Vec(0, 1, 0), (float)(Math.PI/2)));
594                 goal.transform(Matrix.rotate(new Vec(0, 0, 1), (float)(Math.PI/2)));
595                 goal.transform(Matrix.scale(1, 2.2f, 1));
596                 float factor = (float)Math.pow(tile.volume() / goal.volume(), 1.0/3.0);
597                 factor = factor * 0.8f;
598                 goal.transform(Matrix.scale(factor));
599                 fixupGoal(true, false);
600                 //goal.transform(Matrix.translate(new Vec(0.145f, 0, 0)));
601                 fixupGoal(false, true);
602             }});
603             goalMenu.add(new MyMenuItem("Hammerhead Fish") { public void hit() {
604                 loadGoal("fish.stl");
605                 goal.transform(Matrix.rotate(new Vec(0, 1, 0), (float)(Math.PI/2)));
606                 goal.transform(Matrix.rotate(new Vec(0, 0, 1), (float)(3*Math.PI/2)));
607                 float factor = (float)Math.pow(tile.volume() / goal.volume(), 1.0/3.0);
608                 factor *= 0.75f;
609                 goal.transform(Matrix.scale(factor));
610                 fixupGoal();
611             }});
612             goalMenu.add(new MyMenuItem("Vertical Fish") { public void hit() {
613                 loadGoal("fish.stl");
614                 //goal.transform(Matrix.rotate(new Vec(0, 0, 1), (float)(Math.PI/2)));
615                 //goal.transform(Matrix.rotate(new Vec(1, 0, 0), (float)(Math.PI/2)));
616                 float factor = (float)Math.pow(tile.volume() / goal.volume(), 1.0/3.0);
617                 goal.transform(Matrix.scale(factor/1.6f));
618                 fixupGoal();
619             }});
620             goalMenu.add(new MyMenuItem("Torus") { public void hit() {
621                 loadGoal("torus.stl");
622                 float factor = (float)Math.pow(tile.volume() / goal.volume(), 1.0/3.0);
623                 goal.transform(Matrix.scale(factor/2.5f));
624                 goal.transform(Matrix.rotate(new Vec(0, 1, 0), (float)(Math.PI/2)));
625                 fixupGoal();
626             }});
627             tileMenu.add(new MyMenuItem("Hex Brick") { public void hit() {
628                 hexBrick(false, false);
629             }});
630             tileMenu.add(new MyMenuItem("Hex Brick, offset planes") { public void hit() {
631                 hexBrick(true, false);
632             }});
633             tileMenu.add(new MyMenuItem("Hex Brick, rotated") { public void hit() {
634                 hexBrick(false, true);
635             }});
636             tileMenu.add(new MyMenuItem("Temp (do not use)") { public void hit() {
637                 setTile(new Mesh(false));
638                 float width  = (float)0.8;
639                 float depth  = (float)0.08;
640                 float height = (float)0.4;
641
642                 float rshift =   width/2;
643                 float lshift = -(width/2);
644                 float halfup = 0;
645                 //float shift = height/2;
646                 //width = (width*2)/3;
647                 float shift = 0;
648                 transforms = new Matrix[] {
649
650                     Matrix.translate(new Vec(lshift/2,  depth,    -shift)),
651                     Matrix.translate(new Vec(rshift/2,  depth,    -shift)),
652                     Matrix.translate(new Vec(lshift/2, -depth,    -shift)),
653                     Matrix.translate(new Vec(rshift/2, -depth,    -shift)),
654
655                     Matrix.translate(new Vec(lshift,  depth/2,    -shift)),
656                     Matrix.translate(new Vec(rshift,  depth/2,    -shift)),
657                     Matrix.translate(new Vec(lshift, -depth/2,    -shift)),
658                     Matrix.translate(new Vec(rshift, -depth/2,    -shift)),
659
660
661                     /*
662                       Matrix.translate(new Vec(lshift,  depth,    -shift)),
663                       Matrix.translate(new Vec(rshift,  depth,    -shift)),
664                       Matrix.translate(new Vec(lshift, -depth,    -shift)),
665                       Matrix.translate(new Vec(rshift, -depth,    -shift)),
666                     */
667                     /*
668                       Matrix.translate(new Vec(lshift,  depth,    shift)),
669                       Matrix.translate(new Vec(rshift,  depth,    shift)),
670                       Matrix.translate(new Vec(lshift, -depth,    shift)),
671                       Matrix.translate(new Vec(rshift, -depth,    shift)),
672                     */
673                     //Matrix.translate(new Vec(0,  depth,    0)).times(Matrix.rotate(new Vec(0, 0, 1), (float)Math.PI)),
674                     //Matrix.translate(new Vec(0, -depth,    0)).times(Matrix.rotate(new Vec(0, 0, 1), (float)Math.PI)),
675                     //Matrix.translate(new Vec(0,   0,    height)).times(Matrix.rotate(new Vec(0, 0, 1), (float)Math.PI)),
676                     //Matrix.translate(new Vec(0,   0,   -height)).times(Matrix.rotate(new Vec(0, 0, 1), (float)Math.PI)),
677
678                     //Matrix.translate(new Vec(0,  depth, 0)),
679                     //Matrix.translate(new Vec(0, -depth, 0)),
680                     Matrix.translate(new Vec(0,      0,  height)),
681                     Matrix.translate(new Vec(0,      0, -height)),
682
683                     //Matrix.translate(new Vec(lshift,        depth,  height/2)),
684                     //Matrix.translate(new Vec(lshift,        depth, -height/2)),
685                     //Matrix.translate(new Vec(rshift,       -depth,  height/2)),
686                     //Matrix.translate(new Vec(rshift,       -depth, -height/2)),
687                     //Matrix.translate(new Vec(rshift,       0,  height)),
688                     //Matrix.translate(new Vec(rshift,       0, -height)),
689
690                     Matrix.translate(new Vec( width,           0,    0)),
691                     Matrix.translate(new Vec(-width,           0,    0)),
692
693                     Matrix.ONE
694                 };
695                 fixupTile();
696             } });
697             tileMenu.add(new MyMenuItem("Dense Packing (hex)") { public void hit() {
698                 setTile(new Mesh(false));
699                 float width  = (float)3.2;
700                 float depth  = (float)0.32;
701                 float height = (float)1.6;
702                 float unit = 0.4f;
703                 float r = unit/2;
704                 float sin = (float)(unit * Math.sin(Math.PI/3));
705                 float cos = (float)(unit * Math.cos(Math.PI/3));
706                 float x = (float)(r*Math.tan(Math.PI/6));
707                 float z = (float)(r/Math.cos(Math.PI/6));
708                 height = 2*r*(float)Math.sqrt(2f/3f);
709
710                 /*
711                 r *= 0.3f;
712                 cos *= 0.3f;
713                 unit *= 0.3f;
714                 */
715
716                 /*
717                   sin *= 0.3f;
718                   x *= 0.3f;
719                   z *= 0.3f;
720                 */
721                 transforms = new Matrix[] {
722                     Matrix.translate(new Vec(-unit, 0, 0)),
723                     Matrix.translate(new Vec( unit, 0, 0)),
724                     Matrix.translate(new Vec(-cos,  0,  sin)),
725                     Matrix.translate(new Vec( cos,  0,  sin)),
726                     Matrix.translate(new Vec(-cos,  0, -sin)),
727                     Matrix.translate(new Vec( cos,  0, -sin)),
728                     Matrix.translate(new Vec( 0,  height, z)),
729                     Matrix.translate(new Vec(-r,  height, -x)),
730                     Matrix.translate(new Vec( r,  height, -x)),
731                     Matrix.translate(new Vec( 0, -height, -z)),
732                     Matrix.translate(new Vec(-r, -height, x)),
733                     Matrix.translate(new Vec( r, -height, x)),
734                     Matrix.ONE,
735                 };
736                 generateTile(transforms, tile);
737                 fixupTile();
738             } });
739
740             tileMenu.add(new MyMenuItem("Slim Dense Packing (Cubic)") { public void hit() {
741                 setTile(new Mesh(false));
742                 float unit = 0.4f;
743                 float r = unit/2;
744                 float sin = (float)(unit * Math.sin(Math.PI/3));
745                 float cos = (float)(unit * Math.cos(Math.PI/3));
746                 float x = (float)(r*Math.tan(Math.PI/6));
747                 float z = (float)(r/Math.cos(Math.PI/6));
748                 float height = 2*r*(float)Math.sqrt(2f/3f);
749
750                 transforms = new Matrix[] {
751
752                     //Matrix.reflect(new Vec( 0,  height, -z).norm()),
753
754                     Matrix.translate(new Vec(-unit, 0, 0)),
755                     Matrix.translate(new Vec( unit, 0, 0)),
756                     Matrix.translate(new Vec(-cos,  0,  sin)),
757                     Matrix.translate(new Vec( cos,  0,  sin)),
758                     Matrix.translate(new Vec(-cos,  0, -sin)),
759                     Matrix.translate(new Vec( cos,  0, -sin)),
760
761                     Matrix.translate(new Vec( 0,  height, -z)),
762                     Matrix.translate(new Vec(-r,  height,  x)),
763                     Matrix.translate(new Vec( r,  height,  x)),
764                     Matrix.translate(new Vec( 0, -height,  z)),
765                     Matrix.translate(new Vec(-r, -height, -x)),
766                     Matrix.translate(new Vec( r, -height, -x)),
767
768                     /*
769                     Matrix.translate(new Vec( 0,  height, -z)).times(Matrix.rotate(new Vec(0,1,0), (float)Math.PI)),
770                     Matrix.translate(new Vec(-r,  height,  x)).times(Matrix.rotate(new Vec(0,1,0), (float)Math.PI)),
771                     Matrix.translate(new Vec( r,  height,  x)).times(Matrix.rotate(new Vec(0,1,0), (float)Math.PI)),
772                     Matrix.translate(new Vec( 0, -height, -z)).times(Matrix.rotate(new Vec(0,1,0), (float)Math.PI)),
773                     Matrix.translate(new Vec(-r, -height,  x)).times(Matrix.rotate(new Vec(0,1,0), (float)Math.PI)),
774                     Matrix.translate(new Vec( r, -height,  x)).times(Matrix.rotate(new Vec(0,1,0), (float)Math.PI)),
775                     */
776
777                     /*
778                     Matrix.translate(new Vec( 0,  height, -z)).times(Matrix.scale(-1,1,-1)),
779                     Matrix.translate(new Vec(-r,  height,  x)).times(Matrix.scale(-1,1,-1)),
780                     Matrix.translate(new Vec( r,  height,  x)).times(Matrix.scale(-1,1,-1)),
781                     Matrix.translate(new Vec( 0, -height, -z)).times(Matrix.scale(-1,1,-1)),
782                     Matrix.translate(new Vec(-r, -height,  x)).times(Matrix.scale(-1,1,-1)),
783                     Matrix.translate(new Vec( r, -height,  x)).times(Matrix.scale(-1,1,-1)),
784                     */
785                     Matrix.ONE
786                 };
787                 generateTile(transforms, tile);
788
789
790                 transforms = new Matrix[] {
791                     Matrix.reflect(new Vec( 0,  height, -z).norm()),
792                     //Matrix.reflect(new Vec( 0, -height,  z)),
793                     Matrix.ONE,
794                     Matrix.translate(new Vec(-unit, 0, 0)),
795                     Matrix.translate(new Vec( unit, 0, 0)),
796                     /*
797                     Matrix.translate(new Vec(-cos,  0,  sin)),
798                     Matrix.translate(new Vec( cos,  0,  sin)),
799                     Matrix.translate(new Vec(-cos,  0, -sin)),
800                     Matrix.translate(new Vec( cos,  0, -sin)),
801
802                     Matrix.translate(new Vec( 0,  height, -z)),
803                     Matrix.translate(new Vec(-r,  height,  x)),
804                     Matrix.translate(new Vec( r,  height,  x)),
805                     Matrix.translate(new Vec( 0, -height,  z)),
806                     Matrix.translate(new Vec(-r, -height, -x)),
807                     Matrix.translate(new Vec( r, -height, -x)),
808                     */
809                 };
810
811
812                 //for(int i=0; i<transforms.length; i++) transforms[i] = Matrix.translate(m.times(vecs[i]));
813
814                 Matrix m = Matrix.scale(1.9f, 1f ,1);
815                 //Matrix m = Matrix.scale(1f, 2.1f, 1f);
816                 tile.transform(m);
817                 for(int i=0; i<transforms.length; i++)
818                     transforms[i] = preMultiplyTranslationalComponentBy(transforms[i], m);
819
820
821                 fixupTile();
822
823             } });
824             tileMenu.add(new MyMenuItem("Genus-1") { public void hit() {
825                 synchronized(this) {
826                     setTile(new Mesh(false));
827                     Matrix mm = Matrix.scale(0.1f);
828                     float height = 4;
829                     float width  = 4;
830                     float depth  = 1;
831                     // top
832                     quad(tile, mm, 
833                          new Point( 2,  2,  0),
834                          new Point( 1,  1, -1),
835                          new Point(-1,  1, -1),
836                          new Point(-2,  2,  0));
837                     quad(tile, mm, 
838                          new Point(-2,  2,  0),
839                          new Point(-1,  1,  1),
840                          new Point( 1,  1,  1),
841                          new Point( 2,  2,  0));
842                     quad(tile, mm, 
843                          new Point( 1,  1, -1),
844                          new Point( 1,  1,  1),
845                          new Point(-1,  1,  1),
846                          new Point(-1,  1, -1));
847
848                     // bottom
849                     quad(tile, mm, 
850                          new Point(-2, -2,  0),
851                          new Point(-1, -1, -1),
852                          new Point( 1, -1, -1),
853                          new Point( 2, -2,  0));
854                     quad(tile, mm, 
855                          new Point( 2, -2,  0),
856                          new Point( 1, -1,  1),
857                          new Point(-1, -1,  1),
858                          new Point(-2, -2,  0));
859                     quad(tile, mm, 
860                          new Point(-1, -1, -1),
861                          new Point(-1, -1,  1),
862                          new Point( 1, -1,  1),
863                          new Point( 1, -1, -1));
864
865                     // left
866                     quad(tile, mm, 
867                          new Point( 2, -2,  0),
868                          new Point( 1, -1, -1),
869                          new Point( 1,  1, -1),
870                          new Point( 2,  2,  0));
871                     quad(tile, mm, 
872                          new Point( 2,  2,  0),
873                          new Point( 1,  1,  1),
874                          new Point( 1, -1,  1),
875                          new Point( 2, -2,  0));
876                     quad(tile, mm, 
877                          new Point( 1, -1, -1),
878                          new Point( 1, -1,  1),
879                          new Point( 1,  1,  1),
880                          new Point( 1,  1, -1));
881
882                     // bottom
883                     quad(tile, mm, 
884                          new Point(-2,  2,  0),
885                          new Point(-1,  1, -1),
886                          new Point(-1, -1, -1),
887                          new Point(-2, -2,  0));
888                     quad(tile, mm, 
889                          new Point(-2, -2,  0),
890                          new Point(-1, -1,  1),
891                          new Point(-1,  1,  1),
892                          new Point(-2,  2,  0));
893                     quad(tile, mm, 
894                          new Point(-1,  1, -1),
895                          new Point(-1,  1,  1),
896                          new Point(-1, -1,  1),
897                          new Point(-1, -1, -1));
898
899                     height = 4;
900                     width  = 4;
901                     depth  = 1;
902
903
904                     transforms = new Matrix[] {
905                         Matrix.translate(new Vec(0, 0.2f,0))
906                         .times(Matrix.rotate(new Vec(0,1,0), (float)(1*Math.PI/2))),
907
908                         Matrix.translate(new Vec(0,-0.2f,0))
909                         .times(Matrix.rotate(new Vec(0,1,0), (float)(-1*Math.PI/2))),
910
911                         Matrix.translate(new Vec( 0.2f,0,0))
912                         .times(Matrix.rotate(new Vec(1,0,0), (float)(1*Math.PI/2))),
913
914                         Matrix.translate(new Vec(-0.2f,0,0))
915                         .times(Matrix.rotate(new Vec(1,0,0), (float)(-1*Math.PI/2))),
916
917
918                         //Matrix.rotate(new Vec(0,0,1), (float)(1*Math.PI/2)),
919
920                         /*
921                           Matrix.rotate(new Vec(0,0,1), (float)(1*Math.PI/2)),
922
923                           Matrix.rotate(new Vec(0,0,1), (float)(3*Math.PI/2)),
924                           Matrix.rotate(new Vec(1,0,0), (float)(2*Math.PI/2)),
925                         */
926
927                         //Matrix.rotate(new Vec(0,0,1), (float)(2*Math.PI/2)),
928                         //Matrix.scale(1,-1,1),
929
930                         //Matrix.translate(new Vec( 0.2f, 0,0))
931                         //.times(Matrix.rotate(new Vec(0,0,1), (float)( 1*Math.PI/2)))
932                         //.times(Matrix.rotate(new Vec(0,1,0), (float)( 3*Math.PI/2))),
933
934                         //Matrix.translate(new Vec(-0.2f, 0,0))
935                         //.times(Matrix.rotate(new Vec(0,0,1), (float)( 3*Math.PI/2)))
936                         //.times(Matrix.rotate(new Vec(0,1,0), (float)( 3*Math.PI/2))),
937             
938                         //Matrix.rotate(new Vec(0,0,1), (float)( 0*Math.PI/2))
939                         //.times(Matrix.translate(new Vec(0, -0.2f, 0)))
940                         //.times(Matrix.rotate(new Vec(0,1,0), (float)( 1*Math.PI/2))),
941                         //Matrix.rotate(new Vec(0,0,1), (float)( 1*Math.PI/2))
942                         //.times(Matrix.translate(new Vec(0, -0.2f, 0)))
943                         //.times(Matrix.rotate(new Vec(0,1,0), (float)( 1*Math.PI/2))),
944
945                         //Matrix.rotate(new Vec(0,0,1), (float)( 0*Math.PI/2))
946                         //.times(Matrix.translate(new Vec(0, -0.2f, 0)))
947                         //.times(Matrix.rotate(new Vec(0,1,0), (float)( 1*Math.PI/2))),
948                         //Matrix.rotate(new Vec(0,0,1), (float)( 0*Math.PI/2))
949                         //.times(Matrix.translate(new Vec(0, -0.2f, 0)))
950                         //.times(Matrix.rotate(new Vec(0,1,0), (float)( 1*Math.PI/2))),
951
952                         Matrix.ONE,
953                     };
954                     fixupTile();
955                 }}});
956             tileMenu.add(new MyMenuItem("Hammerhead") { public void hit() {
957                 synchronized(this) {
958                     setTile(new Mesh(false));
959                     Matrix mm = Matrix.ONE;
960                     float height1 = .1f;
961                     float height2 = .1f;
962                     float width  = .4f;
963                     float depth  = .1f;
964                     // top
965                     Point a1 = new Point( -width/2,         height2/2, -depth/2);
966                     Point b1 = new Point(        0,         height2/2, -depth/2);
967                     Point c1 = new Point(        0, height1+height2/2, -depth/2);
968                     Point d1 = new Point(  width/2, height1+height2/2, -depth/2);
969                     Point e1 = new Point(  width/2,         height2/2, -depth/2);
970                     Point f1 = new Point(  width/2,        -height2/2, -depth/2);
971                     Point g1 = new Point(  width/2,-height1-height2/2, -depth/2);
972                     Point h1 = new Point(        0,-height1-height2/2, -depth/2);
973                     Point i1 = new Point(        0,        -height2/2, -depth/2);
974                     Point j1 = new Point( -width/2,        -height2/2, -depth/2);
975                     Point a2 = new Point( -width/2,         height2/2,  depth/2);
976                     Point b2 = new Point(        0,         height2/2,  depth/2);
977                     Point c2 = new Point(        0, height1+height2/2,  depth/2);
978                     Point d2 = new Point(  width/2, height1+height2/2,  depth/2);
979                     Point e2 = new Point(  width/2,         height2/2,  depth/2);
980                     Point f2 = new Point(  width/2,        -height2/2,  depth/2);
981                     Point g2 = new Point(  width/2,-height1-height2/2,  depth/2);
982                     Point h2 = new Point(        0,-height1-height2/2,  depth/2);
983                     Point i2 = new Point(        0,        -height2/2,  depth/2);
984                     Point j2 = new Point( -width/2,        -height2/2,  depth/2);
985
986                     quad(tile, mm, a1, b1, i1, j1);
987                     quad(tile, mm, c1, d1, e1, b1);
988                     quad(tile, mm, b1, e1, f1, i1);
989                     quad(tile, mm, i1, f1, g1, h1);
990
991                     quad(tile, mm, j2, i2, b2, a2);
992                     quad(tile, mm, b2, e2, d2, c2);
993                     quad(tile, mm, i2, f2, e2, b2);
994                     quad(tile, mm, h2, g2, f2, i2);
995
996                     quad(tile, mm, d1, d2, e2, e1);
997                     quad(tile, mm, e1, e2, f2, f1);
998                     quad(tile, mm, f1, f2, g2, g1);
999                     quad(tile, mm, h1, g1, g2, h2);
1000                     quad(tile, mm, i2, i1, h1, h2);
1001                     quad(tile, mm, j1, i1, i2, j2);
1002                     quad(tile, mm, a2, a1, j1, j2);
1003                     quad(tile, mm, a1, a2, b2, b1);
1004                     quad(tile, mm, c2, c1, b1, b2);
1005                     quad(tile, mm, c1, c2, d2, d1);
1006
1007                     transforms = new Matrix[] {
1008
1009                         mm.times(Matrix.translate(new Vec(   width,     0, 0))),
1010                         mm.times(Matrix.translate(new Vec(  -width,     0, 0))),
1011                         mm.times(Matrix.translate(new Vec(-width/2, height1+height2, 0))),
1012                         mm.times(Matrix.translate(new Vec( width/2, height1+height2, 0))),
1013                         mm.times(Matrix.translate(new Vec(-width/2,-height1-height2, 0))),
1014                         mm.times(Matrix.translate(new Vec( width/2,-height1-height2, 0))),
1015
1016                         mm.times(Matrix.translate(new Vec(   width/2,     0,  depth))),
1017                         mm.times(Matrix.translate(new Vec(  -width/2,     0,  depth))),
1018                         mm.times(Matrix.translate(new Vec(         0, height1+height2,  depth))),
1019                         mm.times(Matrix.translate(new Vec(         0,-height1-height2,  depth))),
1020
1021                         mm.times(Matrix.translate(new Vec(   width/2,     0, -depth))),
1022                         mm.times(Matrix.translate(new Vec(  -width/2,     0, -depth))),
1023                         mm.times(Matrix.translate(new Vec(         0, height1+height2, -depth))),
1024                         mm.times(Matrix.translate(new Vec(         0,-height1-height2, -depth))),
1025
1026                         Matrix.ONE
1027                     };
1028                     fixupTile();
1029                 }}});
1030             tileMenu.add(new MyMenuItem("Marching Cubes") { public void hit() {
1031                 marchingCubes();
1032             }});
1033             tileMenu.add(new MyMenuItem("Brain") { public void hit() {
1034                 marchingCubes2();
1035             }});
1036
1037             // Finally, add all the menus to the menu bar.
1038             add(tileMenu);
1039             //add(goalMenu);
1040             //add(hideMenu);
1041         }
1042
1043     }
1044
1045     private static Matrix preMultiplyTranslationalComponentBy(Matrix mthis, Matrix mm) {
1046         Vec v = mm.times(mthis.getTranslationalComponent());
1047         return new Matrix(mthis.a, mthis.b, mthis.c, v.x,
1048                           mthis.e, mthis.f, mthis.g, v.y,
1049                           mthis.i, mthis.j, mthis.k, v.z,
1050                           mthis.m, mthis.n, mthis.o, 1);
1051     }
1052
1053 //////////////////////////////////////////////////////////////////////////////
1054     public void breakit() {
1055         int oldverts = verts;
1056         if (verts > 2000 && !force) return;
1057         force = false;
1058         System.out.println("doubling vertices.");
1059         PriorityQueue<Mesh.E> es = new PriorityQueue<Mesh.E>();
1060         for(Mesh.T t : tile) {
1061             es.add(t.e1());
1062             es.add(t.e2());
1063             es.add(t.e3());
1064             Thread.yield();
1065             repaint();
1066         }
1067         for(int i=0; i<Math.min(oldverts,50); i++) {
1068             Mesh.E e = es.poll();
1069             verts++;
1070             e.shatter();
1071             Thread.yield();
1072             repaint();
1073         }
1074         System.out.println("now have " + verts + " vertices; max is 2000");
1075         tile.rebindPoints();
1076     }
1077
1078     public boolean rand(double temp, Mesh.Vertex p) {
1079
1080         p.reComputeErrorAround();
1081         double tile_error = tile.error();
1082         double goal_error = goal.error();
1083
1084         float max = p.averageEdgeLength()/5;
1085         Vec v = new Vec(random.nextFloat(), random.nextFloat(), random.nextFloat());
1086         v = v.norm().times((random.nextFloat() - 0.5f) * max);
1087         Matrix m = Matrix.translate(v);
1088         boolean good = p.move(v, false);
1089         if (!good) { return false; }
1090
1091         double new_tile_error = tile.error();
1092         double new_goal_error = goal.error();
1093         double tile_delta = (new_tile_error - tile_error) / tile_error;
1094         double goal_delta = (new_goal_error - goal_error) / goal_error;
1095         double delta = tile_delta + goal_delta;
1096         double swapProbability = Math.exp((-1 * delta) / temp);
1097
1098         //boolean doSwap = good && (tile_delta <= 0 && goal_delta <= 0);
1099         //boolean doSwap = good && (tile_delta + goal_delta <= 0);
1100         //if (temp < 0.000001) doSwap = good && (tile_delta <= 0 && goal_delta <= 0);
1101         boolean doSwap = good && (Math.random() < swapProbability);
1102
1103         // always move uphill if possible -- potential idea
1104         if (tile_delta <= 0 && goal_delta <= 0 && good) doSwap = true;
1105         if (hillclimb)
1106             doSwap = tile_delta <= 0 && goal_delta <= 0 && good;
1107
1108         if (doSwap) {
1109             tile_error = new_tile_error;
1110             goal_error = new_goal_error;
1111             hits++;
1112             p.goodp = p.p;
1113         } else {
1114             p.move(v.times(-1), true);
1115             misses++;
1116         }
1117         p.reComputeErrorAround();
1118         return true;
1119     }
1120
1121     float hits = 0;
1122     float misses = 0;
1123     public void anneal() throws Exception {
1124         double hightemp = 1;
1125         temp = hightemp;
1126         double last = 10;
1127         boolean seek_upward = false;
1128         double acceptance = 1;
1129         while(true) {
1130             synchronized(this) {
1131                 if (!anneal) { repaint(); Thread.sleep(10); continue; }
1132
1133             double ratio = (hits+misses==0) ? 1 : (hits / (hits+misses));
1134             hits = 0;
1135             misses = 0;
1136             double gamma = 1;
1137             acceptance = (ratio+acceptance)/2;
1138             accepts = (int)(Math.ceil(ratio*100));
1139             temps = (int)(Math.ceil(temp*1000));
1140             vertss = tile.size();
1141             if (breaks > 0) {
1142                 while (breaks>0) {
1143                     breaks--;
1144                     breakit();
1145                 }
1146                 seek_upward = true;
1147                 continue;
1148             } else if (acceptance > 0.96) gamma = 0.1f;
1149             else if (acceptance > 0.9)    gamma = 0.2f;
1150             else if (acceptance > 0.8)    gamma = 0.3f;
1151             else if (acceptance > 0.6)    gamma = 0.4f;
1152             else if (acceptance > 0.3)    gamma = 0.8f;
1153             else if (acceptance > 0.15)   gamma = 0.94f;
1154             else if (acceptance > 0.10)   gamma = 0.98f;
1155
1156             else if (acceptance < 0.01)   breaks++;
1157
1158             if (seek_upward) {
1159                 if (acceptance > 0.25) seek_upward = false;
1160                 else gamma = 2-gamma;
1161             }
1162
1163             if (anneal)
1164                 temp = temp * gamma;
1165
1166
1167             HashSet<Mesh.Vertex> hs = new HashSet<Mesh.Vertex>();
1168             for(Mesh.Vertex p : tile.vertices()) hs.add(p);
1169             Mesh.Vertex[] pts = (Mesh.Vertex[])hs.toArray(new Mesh.Vertex[0]);
1170
1171             int count = 0;
1172             long then = System.currentTimeMillis();
1173             for(int i=0; i<300; i++) {
1174                 if (anneal) {
1175                     Mesh.Vertex v = pts[Math.abs(random.nextInt()) % pts.length];
1176                     if (breaks>0) break;
1177                     if (!rand(temp,v)) { i--; continue; }
1178                     v.recomputeFundamentalQuadricIfStale();
1179                     v.recomputeFundamentalQuadricIfNeighborChanged();
1180                     count++;
1181                 }
1182                 Thread.yield();
1183                 repaint();
1184             }
1185
1186             System.out.println("temp="+temp + " ratio="+(Math.ceil(acceptance*100)) + " " +
1187                                "points_per_second=" +
1188                                (count*1000)/((double)(System.currentTimeMillis()-then)));
1189             for(Mesh.Vertex p : goal.vertices()) {
1190                 p.quadricStale = true;
1191                 p.recomputeFundamentalQuadricIfNeighborChanged();
1192             }
1193             }
1194         }
1195     }
1196
1197
1198     public static void main(String[] s) throws Exception {
1199         JFrame f = new JFrame();
1200         f.setLayout(new BorderLayout());
1201         Main main = new Main(f);
1202         f.add(main, BorderLayout.CENTER);
1203
1204         f.addWindowListener(new WindowListener() {
1205                 public void windowActivated(WindowEvent e) { }
1206                 public void windowClosed(WindowEvent e) { System.exit(0); }
1207                 public void windowClosing(WindowEvent e) { System.exit(0); }
1208                 public void windowDeactivated(WindowEvent e)  { }
1209                 public void windowDeiconified(WindowEvent e)  { }
1210                 public void windowIconified(WindowEvent e)  { }
1211                 public void windowOpened(WindowEvent e) { }
1212             });
1213
1214         f.pack();
1215         f.show();
1216         f.setSize(900, 900);
1217         f.doLayout();
1218
1219         JFrame f2 = new JFrame();
1220         f2.setJMenuBar(main.new MyMenuBar());
1221         f2.setSize(100,100);
1222         f2.show();
1223
1224         main.anneal();
1225     }
1226
1227
1228 }