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