June 9 2011

Glass Tomography: Inhomogeneous Translucent Materials

Goal:

• Reconstruction of the internal structure as well as surface geometry of inhomogeneous translucent materials.

Methods:

• Measure data from visual hull surrounding structure 2 screens and lcd background.
• Acquired data in the form of the magnitude of horizontal and vertical deflections of light rays.
• Very repeatable results 80-90 dB Scan-to-scan SNR.
• To test whether data is accurate representation of the system, could use the data as an environment map and compare it to the actual image of the translucent materials.

Model:

• Solving the equation f( x ) = b where:
• f is a function that represents the image formation model. This may be possible to approximate as a linear model.
• Idea is to match the physical model as close as possible to the real model of the translucent materials.
• x here, the parameter of the model, is a refractive index field.
• There is also reflection, scattering, discretization issues etc., but essentially ignoring these for now to focus on refraction.
• Issue: How do you represent (discretize) the refractive index field x ?
• b is a physical artifact. Measurements of the system. These measurements have noise.
• We're capturing angles of exiting rays from the inhomogeneous translucent material.
• So there's a forward and backward (inverse) problem considering this equation:
• Forward: Given the model, that is, f, what's the output, b ?
• Backward (Inverse): Given the recorded measurement, b , what is the model, that is, f?

• The model parameters:
• Refractive index field. (scalar n)
• Discretization onto a voxel grid. (local basis function)
• Grid resolution. (increasing resolution)
• Regular versus unstructured grid: there may be pay-offs in distributing the grid intelligently versus uniformly. (Performance savings, but tradeoff between x and y coordinates and value of voxels)
• Alternative approach: global wavelet basis function.

• Image Formation Model:
• Ray of geometric operations. ( d/ds(n(d x /ds)) = Gradient(n) )
• To system of 1st order ODE's.
• Solve with Runge-Kutta 4 IVP.
• Adaptive steps. ( d d /ds = Gradient(n))
• ( n(d x /ds) = d )
• Question: Are you using error estimation at each step?
• Answer: Basically do one small step and one large step at each step, compare the error of each, and if the error is comparable between the two, see if we can get away with the smaller step.
• Question: Do you have the ray exit position/direction?
• Answer: We have the ray's exit point from the material. We can then trace to the visual hull. Basically we're only solving for voxels in the visual hull. We're doing 2 plane parameterization for voxelation, that's why we move the lcd display. This gives you the entry and exit position for the ray. It basically gives the same information as light field probes would.
• Note: It takes about 1-2 hours per scan.
• The step size plays a large role in different regions of the material. (Where there are differences in glancing angle)
• Put 20K rays/sec in uniform 1K voxel.
• Question: What is the step size?
• Answer: Currently constructing a regular grid with 0.1-0.3 voxel. Step size is a function of this.
• Question: Do you have any reconstructions?
• Answer: No.
• Question: Do you have a grid?
• Answer: No.
• Question: Did you use a solver?
• Answer: Yes, to reconstruct x (the refractive index field), but not for measurements ( b ).
• Suggestion: You could use MATLAB, compute the Jacobi and vector and use this with the raytracer.

• The equation:
• Find set of parameters that'll minimize the redisual refractive index between 1 and 1.3.
• We expect some kind of Beta/Binary distribution.
• f( x ) is non-linear.
• min || f( x ) - b || (with x belonging to X)
• min || x || such that || f( x ) - b || < sigma (with x belonging to X)
• BPDN.
• Constrained x.
• Cast as optimization problem.
• Trust-Region, etc.
• Implicit linear operation at each iteration of solver.
• Question: Do you calculate the gradient information here?
• Answer: No.
• Suggestion: Well, there's 2 ways of approaching this problem: One, use a general non-linear solver; MATLAB chooses the solver for you automatically. Use function: "fmincon" in MATLAB. But this takes a long time. The second way of solving this problem is if you can compute the Jacobi or evaluation the function, then you can use MATLAB to do this really fast.
• Further Clarification: With "fmincon" (don't specify Jacobi) specify some non-linear function and can explicitly compute Jacobi, or explicitly specify the function that will evaluate the Jacobi at the specific point you want. If you represent the refractive index field as a basis function won't have to do n+1 evaluations to get gradients.
• Suggestion: (James) I have a derivation of ray differentials you can look at.

• Linearization:
• Implicit local approximation in non-linear solver.
• Explicit linearization at each iteration.
• Approach by Cha & Vest. A x = b
• Trace rays, solve, retrace, solve, retrace, ...
• Overlay pyramid resolution iterations or just use AMG.
• If you don't use pyramid will not account for rays curving and passing voxels you didn't expect, this way when you update you'll update the wrong voxels.
• Suggestions: Could you use SPH and adaptively reconstruct regions? Because Essex (Working with Dr. Robert Bridson) may have some work that may help you (method involves moving from calculations on particles to a grid, and back to particles.). It's a mixed particle/grid solver for fluids. It solves on particles, different parts of the system too difficult to solver are solved on a grid, and then it moves back to particles.
• Suggestions: You could start with a low res grid, or start with full system and throw the whole thing into a multigrid solver.
• Suggestions: Would have to be a non-linear multigrid solver, though. If you have an ill-conditioned system the multigrid solver is more robust and pushes the system to smoother representation.
• Suggestions: Non-regular grid is a good idea because your system is sparse.
• Suggestions: Could use regular grid with sparsity constraints.
• Suggestions: Could use regular grid, HOWEVER, there is a real practical problem with storing large systems with a regular grid.
• Suggestions: There are pre-packaged, off the shelf coarse structured grid codes (Samurai, etc.) that are available online. They're cache friendly and still adaptive. Start with coarse grid, and for every voxel need to fit it to a finer grid.
• Suggestions: However, these off the shelf codes are very scene-specific, and involve a lot of fine tuning of parameters to a specific scene setup.

• Advantage of A:
• Fast and easy to compute A x = A^T b
• All machinery of linear algebra will be at our disposal.
• Understand the system properties (Can then analyze the rank and condition number etc.).
• Note that for full resolution system this will be difficult to compute.
• Additional benefit is for proofs.

• Disadvantage of A:
• Sheer size. Back of the envelope calculation: (512^3 * 0.10)^2 will be the number of unknowns. (512^3 * 0.10)^2 * 460 will be the total number of entries. This is roughly 6184752906 entries. If each entry is 4 bytes, then memory required will be about 23.5 Gigabytes. Yikes.
• So memory is an issue.
• May have to go for Matrix-free representation to avoid this.
• For low res system could be okay, but otherwise, probably not.

• Solving for Gradients:
• f(nx) = deltax
• f(ny) = deltay
• f(nz) = deltaz

• Linear Solvers:
• QR (This option we can probably take out, because of the sheer size of the problem.)
• SART
• SPGLI
• Suggestion: Can approximate the hessian using finite difference.
• We want to use LBFGs
• Can you still use it?
• Still need to update the hessian but don't need to compute it.
• Try using MATLAB with non-linear solver full grid with Jacobi in raytracer. See how far you can push this in terms of resolution and size.
• Or could try non-linear solver with outer iterations you can specify. At each step solve the system (but solve for refractive index field, not gradients).
• If you can get gradients, however, can solve non-linear system by parallelizing processes across rays.
• If calculate the accumulated gradient as you trace the ray, most of the computation is done in tracing the ray.

• PSM BRAINSTORMING

• Reconstructing Geophysical Properties:
• For instance, the hot spot under Hawai'i was further away from where they thought it was using this kind of reconstruction.
• Tomographic exploration of underground.
• Methods for generating and recording measurements in this problem space: plane of explosives and measure waves; limited angle tomography (sensors on earth's surface)
• Same principle applies to atmosphere, lakes, etc.
• Are rays still a good model here at large wavelengths?

• Biological Systems:
• Visualization of where nerves/veins are exactly.
• Tomography/Crystallography of biological molecules: investigators are able to take many, many images of the molecule, however, run into a problem in stitching them together to reconstruct the molecule's 3D structure since they're not sure about the direction the image was taken from, the images are noisy, and un-calibrated.

END OF MEETING

-- Main.joelaf - 08 Jun 2011