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     }
 230
 231     public synchronized void breakit() {
 232         if (verts > 800) return;
 233         //while(verts < 800) {
 234         PriorityQueue<Mesh.E> es = new PriorityQueue<Mesh.E>();
 235         for(Mesh.E e : tile.edges()) es.add(e);
 236         for(int i=0; i<8; i++) {
 237             Mesh.E e = es.poll();
 238             verts++;
 239             System.out.println("shatter " + e);
 240             e.shatter();
 241         }
 242         tile.rebindPoints();
 243         //}
 244     }
 245
 246     public synchronized void rand(float temp, Mesh.Vert p) {
 247
 248         //p.reComputeError();
 249         p.reComputeErrorAround();
 250         double tile_score = tile.score();
 251         double goal_score = goal.score();
 252
 253         Vec v = new Vec((random.nextFloat() - (float)0.5) / 200,
 254                         (random.nextFloat() - (float)0.5) / 200,
 255                         (random.nextFloat() - (float)0.5) / 200);
 256         /*
 257         Matrix inv = p.errorQuadric();
 258         v = new Vec(inv.d, inv.h, inv.l).norm().times(1/(float)1000);
 259         */
 260
 261         boolean good = p.move(v);
 262
 263         p.reComputeErrorAround();
 264
 265         double new_tile_score = tile.score();
 266         double new_goal_score = goal.score();
 267         double tile_delta = (new_tile_score - tile_score) / tile_score;
 268         double goal_delta = (new_goal_score - goal_score) / goal_score;
 269         double delta = tile_delta + goal_delta;
 270         double swapProbability = Math.exp((-1 * delta) / temp);
 271         boolean doSwap = good && (Math.random() < swapProbability);
 272         //boolean doSwap = good && (tile_delta <= 0 && goal_delta <= 0);
 273         if (doSwap) {
 274             tile_score = new_tile_score;
 275             goal_score = new_goal_score;
 276             //System.out.println("score: " + tile_score + " / " + goal_score);
 277             hits++;
 278         } else {
 279             p.move(v.times(-1));
 280             misses++;
 281         }
 282     }
 283
 284     float hits = 0;
 285     float misses = 0;
 286     public void anneal() throws Exception {
 287         int verts = 0;
 288         float hightemp = 10;
 289         float temp = hightemp;
 290         float last = 10;
 291         while(true) {
 292             double ratio = (hits+misses==0) ? 1 : (hits / (hits+misses));
 293             System.out.println("temp="+temp + " ratio="+(Math.ceil(ratio*100)));
 294             hits = 0;
 295             misses = 0;
 296             float gamma = 0;
 297             double acceptance = ratio;
 298             if (breaks) {
 299                 breaks = false;
 300                     breakit();
 301                     gamma = 1;
 302                     //gamma = 2;
 303                     //temp = last * 0.8f;
 304                     //last = temp;
 305                     //temp = hightemp;
 306             } else
 307             if (acceptance > 0.96) gamma = 0.6f;
 308             else if (acceptance > 0.9) gamma = 0.7f;
 309             else if (acceptance > 0.8) gamma = 0.75f;
 310             else if (acceptance > 0.6) gamma = 0.8f;
 311             else {
 312                 if (acceptance > 0.3) {
 313                     gamma = 0.9f;
 314                 } else if (acceptance > 0.15) {
 315                     gamma = 0.95f;
 316                 } else {
 317                     breakit();
 318                     //gamma = 2;
 319                     gamma = 1;
 320                     //temp = last * 0.8f;
 321                     //last = temp;
 322                     //temp = hightemp;
 323                 }
 324             }
 325             temp = temp * gamma;
 326
 327             HashSet<Mesh.Vert> hs = new HashSet<Mesh.Vert>();
 328             for(Mesh.Vert p : tile.vertices()) hs.add(p);
 329             for(int i=0; i<10; i++) {
 330                 repaint();
 331                 for(Mesh.Vert v : hs) rand(temp,v);
 332             }
 333             tile.rebuildPointSet();
 334             repaint();
 335             //breakit();
 336             repaint();
 337             goal.unApplyQuadricToNeighborAll();
 338             repaint();
 339             tile.recomputeAllFundamentalQuadrics();
 340             repaint();
 341             goal.applyQuadricToNeighborAll();
 342        }
 343     }
 344
 345     public static void main(String[] s) throws Exception {
 346         StlFile stlf = new StlFile();
 347         stlf.load("simplefish.stl");
 348         Frame f = new Frame();
 349         Main main = new Main(stlf, f);
 350         f.pack();
 351         f.show();
 352         f.setSize(900, 900);
 353         f.doLayout();
 354         main.anneal();
 355     }
 356
 357 }