src/edu/berkeley/qfat/Main.java

   1 package edu.berkeley.qfat;
   2 import java.awt.*;
   3 import java.awt.event.*;
   4 import javax.swing.*;
   5 import javax.media.opengl.*;
   6 import javax.media.opengl.glu.*;
   7 import java.util.*;
   8 import edu.berkeley.qfat.geom.*;
   9 import edu.berkeley.qfat.geom.Point;
  10
  11 // TO DO:
  12 // - real anneal
  13 // - solve self-intersection problem
  14 // - get a better test model?
  15
  16 // FIXME: re-orient goal (how?)
  17
  18 public class Main extends MeshViewer {
  19
  20     public static int verts = 0;
  21
  22     public static final Random random = new Random();
  23
  24     /** magnification factor */
  25     private static final float MAG = 1;
  26
  27     public Main(StlFile stlf, Frame f) {
  28         super(f);
  29
  30         for(int i=0; i<stlf.coordArray.length; i+=3) {
  31             Point p0 = new Point(stlf.coordArray[i+0].x * MAG, stlf.coordArray[i+0].y * MAG, stlf.coordArray[i+0].z * MAG);
  32             Point p1 = new Point(stlf.coordArray[i+1].x * MAG, stlf.coordArray[i+1].y * MAG, stlf.coordArray[i+1].z * MAG);
  33             Point p2 = new Point(stlf.coordArray[i+2].x * MAG, stlf.coordArray[i+2].y * MAG, stlf.coordArray[i+2].z * MAG);
  34             Vec n  = new Vec(stlf.normArray[i/3].x * MAG, stlf.normArray[i/3].y  * MAG, stlf.normArray[i/3].z * MAG);
  35             Mesh.T t  = goal.newT(p0, p1, p2, n, 0);
  36         }
  37
  38         // rotate to align major axis -- this probably needs to be done by a human.
  39         goal.transform(new Matrix(new Vec(0, 0, 1), (float)(Math.PI/2)));
  40
  41         float goal_width  = goal.diagonal().dot(new Vec(1, 0, 0));
  42         float goal_height = goal.diagonal().dot(new Vec(0, 1, 0));
  43         float goal_depth  = goal.diagonal().dot(new Vec(0, 0, 1));
  44
  45         /*
  46         float width  = (float)0.6;
  47         float height = (float)0.08;
  48         float depth  = (float)0.3;
  49         */
  50         float width  = (float)0.6;
  51         float depth  = (float)0.08;
  52         float height = (float)0.3;
  53
  54         float rshift = width/2;
  55         float lshift = -(width/2);
  56
  57         //float halfup = height/2;
  58         float halfup = 0;
  59
  60         translations = new Matrix[] {
  61
  62             new Matrix(new Vec(lshift,  depth,    halfup)),
  63             new Matrix(new Vec(rshift,  depth,    halfup)),
  64             new Matrix(new Vec(lshift, -depth,    halfup)),
  65             new Matrix(new Vec(rshift, -depth,    halfup)),
  66
  67             new Matrix(new Vec(lshift,       0,  height)),
  68             new Matrix(new Vec(rshift,       0,  height)),
  69             new Matrix(new Vec(lshift,       0, -height)),
  70             new Matrix(new Vec(rshift,       0, -height)),
  71
  72             new Matrix(new Vec( width,           0,    0)),
  73             new Matrix(new Vec(-width,           0,    0)),
  74
  75         };
  76
  77
  78         Point ltf = new Point(lshift,  (depth/2),  (height/2));
  79         Point mtf = new Point( 0.0,    (depth/2),  (height/2));
  80         Point rtf = new Point(rshift,  (depth/2),  (height/2));
  81         Point lbf = new Point(lshift, -(depth/2),  (height/2));
  82         Point mbf = new Point( 0.0,   -(depth/2),  (height/2));
  83         Point rbf = new Point(rshift, -(depth/2),  (height/2));
  84
  85         Point ltc = new Point(lshift,  (depth/2), 0);
  86         Point mtc = new Point( 0.0,    (depth/2), 0);
  87         Point rtc = new Point(rshift,  (depth/2), 0);
  88         Point lbc = new Point(lshift, -(depth/2), 0);
  89         Point mbc = new Point( 0.0,   -(depth/2), 0);
  90         Point rbc = new Point(rshift, -(depth/2), 0);
  91
  92         Point ltn = new Point(lshift,  (depth/2), -(height/2));
  93         Point mtn = new Point( 0.0,    (depth/2), -(height/2));
  94         Point rtn = new Point(rshift,  (depth/2), -(height/2));
  95         Point lbn = new Point(lshift, -(depth/2), -(height/2));
  96         Point mbn = new Point( 0.0,   -(depth/2), -(height/2));
  97         Point rbn = new Point(rshift, -(depth/2), -(height/2));
  98
  99
 100         Point[] points = new Point[] {
 101             ltf,
 102             mtf,
 103             rtf,
 104             lbf,
 105             mbf,
 106             rbf,
 107
 108             ltc,
 109             mtc,
 110             rtc,
 111             lbc,
 112             mbc,
 113             rbc,
 114
 115             ltn,
 116             mtn,
 117             rtn,
 118             lbn,
 119             mbn,
 120             rbn
 121         };
 122
 123
 124         // top
 125         tile.newT(ltf, mtf, mtc, null, 1);
 126         tile.newT(mtc, ltc, ltf, null, 1);
 127         tile.newT(mtf, rtf, rtc, null, 1);
 128         tile.newT(rtc, mtc, mtf, null, 1);
 129
 130         tile.newT(ltc, mtc, mtn, null, 1);
 131         tile.newT(mtn, ltn, ltc, null, 1);
 132         tile.newT(mtc, rtc, rtn, null, 1);
 133         tile.newT(rtn, mtn, mtc, null, 1);
 134
 135         // bottom (swap normals)
 136         tile.newT(mbf, lbf, mbc, null, 2);
 137         tile.newT(lbc, mbc, lbf, null, 2);
 138         tile.newT(rbf, mbf, rbc, null, 2);
 139         tile.newT(mbc, rbc, mbf, null, 2);
 140
 141         tile.newT(mbc, lbc, mbn, null, 2);
 142         tile.newT(lbn, mbn, lbc, null, 2);
 143
 144         tile.newT(rbc, mbc, rbn, null, 2);
 145         tile.newT(mbn, rbn, mbc, null, 2);
 146
 147
 148         // left
 149         tile.newT(ltf, ltc, lbc, null, 3);
 150         tile.newT(lbc, lbf, ltf, null, 3);
 151         tile.newT(ltc, ltn, lbn, null, 3);
 152         tile.newT(lbn, lbc, ltc, null, 3);
 153
 154         // right (swap normals)
 155         tile.newT(rtc, rtf, rbc, null, 4);
 156         tile.newT(rbf, rbc, rtf, null, 4);
 157         tile.newT(rtn, rtc, rbn, null, 4);
 158         tile.newT(rbc, rbn, rtc, null, 4);
 159
 160         // front
 161         tile.newT(ltn, mtn, mbn, null, 5);
 162         tile.newT(ltn, mbn, lbn, null, 5);
 163         tile.newT(mtn, rtn, rbn, null, 5);
 164         tile.newT(mtn, rbn, mbn, null, 5);
 165
 166         // back
 167         tile.newT(mtf, ltf, mbf, null, 6);
 168         tile.newT(mbf, ltf, lbf, null, 6);
 169         tile.newT(rtf, mtf, rbf, null, 6);
 170         tile.newT(rbf, mtf, mbf, null, 6);
 171
 172         for(Matrix m : translations) {
 173             for(Mesh.T t1 : tile) {
 174                 for(Mesh.T t2 : tile) {
 175                     if (t1==t2) continue;
 176
 177                     if ((t1.v1().p.times(m).minus(t2.v1().p).mag() < Mesh.EPSILON) &&
 178                         (t1.v2().p.times(m).minus(t2.v3().p).mag() < Mesh.EPSILON) &&
 179                         (t1.v3().p.times(m).minus(t2.v2().p).mag() < Mesh.EPSILON)) {
 180                         t1.e1().bindEdge(t2.e3());
 181                         t1.e2().bindEdge(t2.e2());
 182                         t1.e3().bindEdge(t2.e1());
 183                     }
 184                     if ((t1.v2().p.times(m).minus(t2.v1().p).mag() < Mesh.EPSILON) &&
 185                         (t1.v3().p.times(m).minus(t2.v3().p).mag() < Mesh.EPSILON) &&
 186                         (t1.v1().p.times(m).minus(t2.v2().p).mag() < Mesh.EPSILON)) {
 187                         t1.e2().bindEdge(t2.e3());
 188                         t1.e3().bindEdge(t2.e2());
 189                         t1.e1().bindEdge(t2.e1());
 190                     }
 191                     if ((t1.v3().p.times(m).minus(t2.v1().p).mag() < Mesh.EPSILON) &&
 192                         (t1.v1().p.times(m).minus(t2.v3().p).mag() < Mesh.EPSILON) &&
 193                         (t1.v2().p.times(m).minus(t2.v2().p).mag() < Mesh.EPSILON)) {
 194                         t1.e3().bindEdge(t2.e3());
 195                         t1.e1().bindEdge(t2.e2());
 196                         t1.e2().bindEdge(t2.e1());
 197                     }
 198
 199                 }
 200             }
 201         }
 202
 203         //xMesh.Vert mid = lbf.getE(mbn).shatter();
 204
 205         // rescale to match volume
 206         float factor = (float)Math.pow(tile.volume() / goal.volume(), 1.0/3.0);
 207         goal.transform(new Matrix(factor));
 208
 209         // translate to match centroid
 210         goal.transform(new Matrix(tile.centroid().minus(goal.centroid())));
 211
 212         //tx.e2.shatter();
 213         //tx.e3.shatter();
 214
 215
 216         tile.rebindPoints();
 217
 218         //mid.move(new Vec((float)0,0,(float)-0.05));
 219         //ltn.move(new Vec((float)0,0,(float)-0.05));
 220
 221         //mtf.move(new Vec(0, (float)-0.05, (float)0.05));
 222
 223
 224         System.out.println("tile volume: " + tile.volume());
 225         System.out.println("goal volume: " + goal.volume());
 226
 227         tile.score_against = goal;
 228         goal.score_against = tile;
 229         tile.tilemesh = true;
 230     }
 231
 232     public synchronized void breakit() {
 233         if (verts > 800) return;
 234         //while(verts < 800) {
 235         PriorityQueue<Mesh.E> es = new PriorityQueue<Mesh.E>();
 236         for(Mesh.E e : tile.edges()) es.add(e);
 237         for(int i=0; i<8; i++) {
 238             Mesh.E e = es.poll();
 239             verts++;
 240             System.out.println("shatter " + e);
 241             e.shatter();
 242         }
 243         tile.rebindPoints();
 244         //}
 245     }
 246
 247     public synchronized void rand(float temp, Mesh.Vert p) {
 248
 249         //p.reComputeError();
 250         p.reComputeErrorAround();
 251         double tile_score = tile.score();
 252         double goal_score = goal.score();
 253
 254         Vec v = new Vec((random.nextFloat() - (float)0.5) / 200,
 255                         (random.nextFloat() - (float)0.5) / 200,
 256                         (random.nextFloat() - (float)0.5) / 200);
 257         /*
 258         Matrix inv = p.errorQuadric();
 259         v = new Vec(inv.d, inv.h, inv.l).norm().times(1/(float)1000);
 260         */
 261
 262         boolean good = p.move(v);
 263
 264         p.reComputeErrorAround();
 265
 266         double new_tile_score = tile.score();
 267         double new_goal_score = goal.score();
 268         double tile_delta = (new_tile_score - tile_score) / tile_score;
 269         double goal_delta = (new_goal_score - goal_score) / goal_score;
 270         double delta = tile_delta + goal_delta;
 271         double swapProbability = Math.exp((-1 * delta) / temp);
 272         boolean doSwap = good && (Math.random() < swapProbability);
 273         //boolean doSwap = good && (tile_delta <= 0 && goal_delta <= 0);
 274         if (doSwap) {
 275             tile_score = new_tile_score;
 276             goal_score = new_goal_score;
 277             //System.out.println("score: " + tile_score + " / " + goal_score);
 278             hits++;
 279         } else {
 280             p.move(v.times(-1));
 281             misses++;
 282         }
 283     }
 284
 285     float hits = 0;
 286     float misses = 0;
 287     public void anneal() throws Exception {
 288         int verts = 0;
 289         float hightemp = 10;
 290         float temp = hightemp;
 291         float last = 10;
 292         while(true) {
 293             double ratio = (hits+misses==0) ? 1 : (hits / (hits+misses));
 294             System.out.println("temp="+temp + " ratio="+(Math.ceil(ratio*100)));
 295             hits = 0;
 296             misses = 0;
 297             float gamma = 0;
 298             double acceptance = ratio;
 299             if (breaks) {
 300                 breaks = false;
 301                     breakit();
 302                     //gamma = 1;
 303                     gamma = 2;
 304                     //temp = last * 0.8f;
 305                     //last = temp;
 306                     //temp = hightemp;
 307             } else
 308             if (acceptance > 0.96) gamma = 0.6f;
 309             else if (acceptance > 0.9) gamma = 0.7f;
 310             else if (acceptance > 0.8) gamma = 0.75f;
 311             else if (acceptance > 0.6) gamma = 0.8f;
 312             else {
 313                 if (acceptance > 0.3) {
 314                     gamma = 0.9f;
 315                 } else if (acceptance > 0.15) {
 316                     gamma = 0.95f;
 317                 } else {
 318                     breakit();
 319                     gamma = 2;
 320                     //gamma = 1;
 321                     //temp = last * 0.8f;
 322                     //last = temp;
 323                     //temp = hightemp;
 324                 }
 325             }
 326             temp = temp * gamma;
 327
 328             HashSet<Mesh.Vert> hs = new HashSet<Mesh.Vert>();
 329             for(Mesh.Vert p : tile.vertices()) hs.add(p);
 330             for(int i=0; i<10; i++) {
 331                 repaint();
 332                 for(Mesh.Vert v : hs) rand(temp,v);
 333             }
 334             tile.rebuildPointSet();
 335             repaint();
 336             //breakit();
 337             repaint();
 338             goal.unApplyQuadricToNeighborAll();
 339             repaint();
 340             tile.recomputeAllFundamentalQuadrics();
 341             repaint();
 342             goal.applyQuadricToNeighborAll();
 343        }
 344     }
 345
 346     public static void main(String[] s) throws Exception {
 347         StlFile stlf = new StlFile();
 348         stlf.load("simplefish.stl");
 349         Frame f = new Frame();
 350         Main main = new Main(stlf, f);
 351         f.pack();
 352         f.show();
 353         f.setSize(900, 900);
 354         f.doLayout();
 355         main.anneal();
 356     }
 357
 358 }