Add java code of DistortionMeshProducer for lens usage, as included in Oculus VR RIFT SDK 0.3.2 (w/ dependency on JOGL)

This code has been added to the oculusvr-sdk sub-module
due to the unclear license situation.

While including these parts w/ the whole SDK,
we do not violate the RIFT SDK licence - to our knowledge.

+++

Code implements JOGL's
'jogamp.opengl.util.stereo.DistortionMesh.Producer' interface.

The implementation heavily uses RIFT SDK mathematical equations,
comments, algorithms and constants (DK1).

RIFT SDK C++ code has been partially ported to Java.

Purpose it to add a development testbed for JOGL
without requiring to bind to the native SDL library
and to experiment w/ some parameters.

This module may be used for educational purposes with JOGL's
  'jogamp.opengl.util.stereo.GenericStereoDevice'

If 'compiled in' with the SDK, reflection it utilized
in one of GenericStereoDevice configuration to load this functionality.

This module is not useful for real VR integration,
since it missies proper parameter detection
as well as sensor support.

For real VR integration, the fully fledged
'jogamp.opengl.oculusvr.OVRStereoDevice' shall be used instead.

The latter is part of JOGL and binds dynamically into the RIFT SDK
native library, which is provided within this sub-module (oculusvr-sdk).

3 files changed, 1136 insertions, 0 deletions

diff --git a/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/DistortionMeshProducer.java b/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/DistortionMeshProducer.java
new file mode 100644
index 0000000..aa64f6c
--- /dev/null
+++ b/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/DistortionMeshProducer.java
@@ -0,0 +1,358 @@
+/**
+ * Copyright 2014 JogAmp Community. All rights reserved.
+ *
+ * Redistribution and use in source and binary forms, with or without modification, are
+ * permitted provided that the following conditions are met:
+ *
+ *    1. Redistributions of source code must retain the above copyright notice, this list of
+ *       conditions and the following disclaimer.
+ *
+ *    2. Redistributions in binary form must reproduce the above copyright notice, this list
+ *       of conditions and the following disclaimer in the documentation and/or other materials
+ *       provided with the distribution.
+ *
+ *    3. Compliance with Oculus VR RIFT SDK LICENSE (see below)
+ *
+ * THIS SOFTWARE IS PROVIDED BY JogAmp Community ``AS IS'' AND ANY EXPRESS OR IMPLIED
+ * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
+ * FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL JogAmp Community OR
+ * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
+ * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+ * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
+ * ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
+ * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ *
+ * The views and conclusions contained in the software and documentation are those of the
+ * authors and should not be interpreted as representing official policies, either expressed
+ * or implied, of JogAmp Community.
+ *
+ * XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
+ *
+ * This file contains mathematical equations, comments and algorithms
+ * used in the Oculus VR RIFT SDK 0.3.2.
+ *
+ * Due to unknown legal status, the 'requested' Copyright tag and disclaimer
+ * below has been added.
+ *
+ * XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
+ *
+ * Copyright © 2014 Oculus VR, Inc. All rights reserved
+ *
+ * Oculus VR, Inc. Software Development Kit License Agreement
+ *
+ * Human-Readable Summary*:
+ *
+ * You are Free to:
+ *
+ * Use, modify, and distribute the Oculus VR Rift SDK in source and binary
+ * form with your applications/software.
+ *
+ * With the Following Restrictions:
+ *
+ * You can only distribute or re-distribute the source code to LibOVR in
+ * whole, not in part.
+ *
+ * Modifications to the Oculus VR Rift SDK in source or binary form must
+ * be shared with Oculus VR.
+ *
+ * If your applications cause health and safety issues, you may lose your
+ * right to use the Oculus VR Rift SDK, including LibOVR.
+ *
+ * The Oculus VR Rift SDK may not be used to interface with unapproved commercial
+ * virtual reality mobile or non-mobile products or hardware.
+
+ * * - This human-readable Summary is not a license. It is simply a convenient
+ * reference for understanding the full Oculus VR Rift SDK License Agreement.
+ * The Summary is written as a user-friendly interface to the full Oculus VR Rift
+ * SDK License below. This Summary itself has no legal value, and its contents do
+ * not appear in the actual license.
+ *
+ * Full-length Legal Copy may be found at:
+ *   http://www.oculusvr.com/licenses/LICENSE-3.1
+ *   http://jogamp.org/git/?p=oculusvr-sdk.git;a=blob;f=LICENSE.txt;hb=HEAD
+ *   Or within this repository: oculusvr-sdk/LICENSE.txt
+ *
+ * THIS RIFT SDK AND ANY COMPONENT THEREOF IS PROVIDED BY OCULUS VR AND
+ * ITS CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING,
+ * BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
+ * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL OCULUS VR AS THE
+ * COPYRIGHT OWNER OR ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
+ * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+ * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,
+ * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+ * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS RIFT
+ * SDK OR THE RIFT SDK DERIVATIVES, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ */
+package jogamp.opengl.oculusvr.stereo.lense;
+
+import jogamp.opengl.util.stereo.DistortionMesh;
+import jogamp.opengl.util.stereo.GenericStereoDevice;
+import jogamp.opengl.util.stereo.ScaleAndOffset2D;
+import jogamp.opengl.util.stereo.DistortionMesh.DistortionVertex;
+
+import com.jogamp.opengl.math.FloatUtil;
+import com.jogamp.opengl.math.VectorUtil;
+import com.jogamp.opengl.util.stereo.EyeParameter;
+
+public class DistortionMeshProducer implements DistortionMesh.Producer {
+    private DistortionSpec[] distortionSpecs;
+    private GenericStereoDevice.Config deviceConfig;
+    private final float[] eyeReliefInMeters;
+
+    public DistortionMeshProducer() {
+        this.distortionSpecs = null;
+        this.deviceConfig = null;
+        this.eyeReliefInMeters = new float[] { 0, 0 };
+    }
+
+    @Override
+    public void init(final GenericStereoDevice.Config deviceConfig, final float[] eyeReliefInMeters) {
+        if( this.deviceConfig != deviceConfig ||
+            this.eyeReliefInMeters[0] != eyeReliefInMeters[0] ||
+            this.eyeReliefInMeters[1] != eyeReliefInMeters[1] ) {
+            System.arraycopy(eyeReliefInMeters, 0, this.eyeReliefInMeters, 0, 2);
+            this.distortionSpecs = DistortionSpec.CalculateDistortionSpec(deviceConfig, eyeReliefInMeters);
+            this.deviceConfig = deviceConfig;
+        }
+    }
+
+    // Pow2 for the Morton order to work!
+    // 4 is too low - it is easy to see the 'wobbles' in the HMD.
+    // 5 is close but you can see pixel differences with even/odd frame checking.
+    // 6 is indistinguishable on a monitor on even/odd frames.
+    private static final int DMA_GridSizeLog2   = 6;
+    private static final int DMA_GridSize       = 1<<DMA_GridSizeLog2;
+    private static final int DMA_NumVertsPerEye = (DMA_GridSize+1)*(DMA_GridSize+1);
+    private static final int DMA_NumTrisPerEye  = (DMA_GridSize)*(DMA_GridSize)*2;
+
+    private static float[] TransformTanFovSpaceToScreenNDC(final DistortionSpec distortion,
+                                                           final float[] tanEyeAngle, final boolean usePolyApprox /*= false*/ ) {
+        final float tanEyeAngleRadius = VectorUtil.normVec2(tanEyeAngle);
+        float tanEyeAngleDistortedRadius = distortion.lens.DistortionFnInverseApprox ( tanEyeAngleRadius );
+        if ( !usePolyApprox ) {
+            tanEyeAngleDistortedRadius = distortion.lens.DistortionFnInverse ( tanEyeAngleRadius );
+        }
+        final float[] vec2Tmp1 = new float[2];
+        final float[] tanEyeAngleDistorted;
+        if ( tanEyeAngleRadius > 0.0f ) {
+            tanEyeAngleDistorted = VectorUtil.scaleVec2(vec2Tmp1, tanEyeAngle, tanEyeAngleDistortedRadius / tanEyeAngleRadius );
+        } else {
+            tanEyeAngleDistorted = tanEyeAngle;
+        }
+
+        final float[] framebufferNDC =
+                VectorUtil.addVec2(vec2Tmp1,
+                        VectorUtil.divVec2(vec2Tmp1, tanEyeAngleDistorted, distortion.tanEyeAngleScale),
+                        distortion.lensCenter );
+        return framebufferNDC;
+    }
+
+    private static void TransformScreenNDCToTanFovSpaceChroma (final float[] resultR, final float[] resultG, final float[] resultB,
+                                                               final DistortionSpec distortion, final float[] framebufferNDC) {
+        // Scale to TanHalfFov space, but still distorted.
+        final float[] vec2Tmp1 = new float[2];
+        final float[] tanEyeAngleDistorted =
+                VectorUtil.scaleVec2(vec2Tmp1, VectorUtil.subVec2(vec2Tmp1, framebufferNDC, distortion.lensCenter), distortion.tanEyeAngleScale);
+        // Distort.
+        final float radiusSquared = ( tanEyeAngleDistorted[0] * tanEyeAngleDistorted[0] )
+                            + ( tanEyeAngleDistorted[1] * tanEyeAngleDistorted[1] );
+        final float[] distortionScales = distortion.lens.DistortionFnScaleRadiusSquaredChroma (radiusSquared);
+        VectorUtil.scaleVec2(resultR, tanEyeAngleDistorted, distortionScales[0]);
+        VectorUtil.scaleVec2(resultG, tanEyeAngleDistorted, distortionScales[1]);
+        VectorUtil.scaleVec2(resultB, tanEyeAngleDistorted, distortionScales[2]);
+    }
+
+    public final DistortionMesh create(final EyeParameter eyeParam, final int distortionBits) {
+        // Find the mapping from TanAngle space to target NDC space.
+        final ScaleAndOffset2D  eyeToSourceNDC = new ScaleAndOffset2D(eyeParam.fovhv);
+        final boolean rightEye = 1 == eyeParam.number;
+
+        // When does the fade-to-black edge start? Chosen heuristically.
+        final float fadeOutBorderFraction = 0.075f;
+
+
+        // Populate vertex buffer info
+        float xOffset = 0.0f;
+        @SuppressWarnings("unused")
+        float uOffset = 0.0f;
+
+        if (rightEye)
+        {
+            xOffset = 1.0f;
+            uOffset = 0.5f;
+        }
+
+        // First pass - build up raw vertex data.
+        final DistortionVertex[] vertices =  new DistortionVertex[DMA_NumVertsPerEye];
+        int currentVertex = 0;
+        final float[] tanEyeAnglesR = new float[2], tanEyeAnglesG = new float[2], tanEyeAnglesB = new float[2];
+
+        for ( int y = 0; y <= DMA_GridSize; y++ ) {
+            for ( int x = 0; x <= DMA_GridSize; x++ ) {
+                final float[] pcurVert = new float[DistortionVertex.def_total_size];
+                int idx = 0;
+
+                final float[] sourceCoordNDC = new float[2];
+                // NDC texture coords [-1,+1]
+                sourceCoordNDC[0] = 2.0f * ( (float)x / (float)DMA_GridSize ) - 1.0f;
+                sourceCoordNDC[1] = 2.0f * ( (float)y / (float)DMA_GridSize ) - 1.0f;
+                final float[] tanEyeAngle = eyeToSourceNDC.convertToTanFovSpace(sourceCoordNDC);
+
+                // This is the function that does the really heavy lifting.
+                final float[] screenNDC = TransformTanFovSpaceToScreenNDC ( distortionSpecs[eyeParam.number], tanEyeAngle, false );
+
+                // VIGNETTE - Fade out at texture edges.
+                float edgeFadeIn = ( 1.0f / fadeOutBorderFraction ) *
+                                   ( 1.0f - Math.max( FloatUtil.abs( sourceCoordNDC[0] ), FloatUtil.abs( sourceCoordNDC[0] ) ) );
+                // Also fade out at screen edges.
+                final float edgeFadeInScreen = ( 2.0f / fadeOutBorderFraction ) *
+                                               ( 1.0f - Math.max( FloatUtil.abs( screenNDC[0] ), FloatUtil.abs( screenNDC[1] ) ) );
+                edgeFadeIn = Math.min( edgeFadeInScreen, edgeFadeIn );
+
+                // TIMEWARP
+                float timewarpLerp;
+                {
+                    switch ( deviceConfig.shutterType )
+                    {
+                    case Global:
+                        timewarpLerp = 0.0f;
+                        break;
+                    case RollingLeftToRight:
+                        // Retrace is left to right - left eye goes 0.0 -> 0.5, then right goes 0.5 -> 1.0
+                        timewarpLerp = screenNDC[0] * 0.25f + 0.25f;
+                        if (rightEye)
+                        {
+                            timewarpLerp += 0.5f;
+                        }
+                        break;
+                    case RollingRightToLeft:
+                        // Retrace is right to left - right eye goes 0.0 -> 0.5, then left goes 0.5 -> 1.0
+                        timewarpLerp = 0.75f - screenNDC[0] * 0.25f;
+                        if (rightEye)
+                        {
+                            timewarpLerp -= 0.5f;
+                        }
+                        break;
+                    case RollingTopToBottom:
+                        // Retrace is top to bottom on both eyes at the same time.
+                        timewarpLerp = screenNDC[1] * 0.5f + 0.5f;
+                        break;
+                    default:
+                        timewarpLerp = 0.0f;
+                        break;
+                    }
+
+                }
+
+                // We then need RGB UVs. Since chromatic aberration is generated from screen coords, not
+                // directly from texture NDCs, we can't just use tanEyeAngle, we need to go the long way round.
+                TransformScreenNDCToTanFovSpaceChroma ( tanEyeAnglesR, tanEyeAnglesG, tanEyeAnglesB,
+                                                        distortionSpecs[eyeParam.number], screenNDC );
+
+                // Don't let verts overlap to the other eye.
+                screenNDC[0] = Math.max ( -1.0f, Math.min( screenNDC[0], 1.0f ) );
+                screenNDC[1] = Math.max ( -1.0f, Math.min( screenNDC[1], 1.0f ) );
+
+                // POS
+                pcurVert[idx++] = 0.5f * screenNDC[0] - 0.5f + xOffset;
+                pcurVert[idx++] = -screenNDC[1];
+
+                // VIGNETTE
+                pcurVert[idx++] = Math.max( 0.0f, Math.min( edgeFadeIn, 1.0f ) );
+
+                // TIMEWARP
+                pcurVert[idx++] = timewarpLerp;
+
+                // Chroma RGB
+                pcurVert[idx++] = tanEyeAnglesR[0];
+                pcurVert[idx++] = tanEyeAnglesR[1];
+                pcurVert[idx++] = tanEyeAnglesG[0];
+                pcurVert[idx++] = tanEyeAnglesG[1];
+                pcurVert[idx++] = tanEyeAnglesB[0];
+                pcurVert[idx++] = tanEyeAnglesB[1];
+
+                vertices[currentVertex++] = new DistortionVertex(
+                        pcurVert, DistortionVertex.def_pos_size,
+                        DistortionVertex.def_vignetteFactor_size, DistortionVertex.def_timewarpFactor_size, DistortionVertex.def_texR_size,
+                        DistortionVertex.def_texG_size, DistortionVertex.def_texB_size);
+            }
+        }
+
+
+        // Populate index buffer info
+        final short[] indices = new short[DMA_NumTrisPerEye*3];
+        int idx = 0;
+
+        for ( int triNum = 0; triNum < DMA_GridSize * DMA_GridSize; triNum++ )
+        {
+            // Use a Morton order to help locality of FB, texture and vertex cache.
+            // (0.325ms raster order -> 0.257ms Morton order)
+            final int x = ( ( triNum & 0x0001 ) >> 0 ) |
+                          ( ( triNum & 0x0004 ) >> 1 ) |
+                          ( ( triNum & 0x0010 ) >> 2 ) |
+                          ( ( triNum & 0x0040 ) >> 3 ) |
+                          ( ( triNum & 0x0100 ) >> 4 ) |
+                          ( ( triNum & 0x0400 ) >> 5 ) |
+                          ( ( triNum & 0x1000 ) >> 6 ) |
+                          ( ( triNum & 0x4000 ) >> 7 );
+
+            final int y = ( ( triNum & 0x0002 ) >> 1 ) |
+                          ( ( triNum & 0x0008 ) >> 2 ) |
+                          ( ( triNum & 0x0020 ) >> 3 ) |
+                          ( ( triNum & 0x0080 ) >> 4 ) |
+                          ( ( triNum & 0x0200 ) >> 5 ) |
+                          ( ( triNum & 0x0800 ) >> 6 ) |
+                          ( ( triNum & 0x2000 ) >> 7 ) |
+                          ( ( triNum & 0x8000 ) >> 8 );
+
+            final int FirstVertex = x * (DMA_GridSize+1) + y;
+            // Another twist - we want the top-left and bottom-right quadrants to
+            // have the triangles split one way, the other two split the other.
+            // +---+---+---+---+
+            // |  /|  /|\  |\  |
+            // | / | / | \ | \ |
+            // |/  |/  |  \|  \|
+            // +---+---+---+---+
+            // |  /|  /|\  |\  |
+            // | / | / | \ | \ |
+            // |/  |/  |  \|  \|
+            // +---+---+---+---+
+            // |\  |\  |  /|  /|
+            // | \ | \ | / | / |
+            // |  \|  \|/  |/  |
+            // +---+---+---+---+
+            // |\  |\  |  /|  /|
+            // | \ | \ | / | / |
+            // |  \|  \|/  |/  |
+            // +---+---+---+---+
+            // This way triangle edges don't span long distances over the distortion function,
+            // so linear interpolation works better & we can use fewer tris.
+            if ( ( x < DMA_GridSize/2 ) != ( y < DMA_GridSize/2 ) )       // != is logical XOR
+            {
+                indices[idx++] = (short)(FirstVertex);
+                indices[idx++] = (short)(FirstVertex+1);
+                indices[idx++] = (short)(FirstVertex+(DMA_GridSize+1)+1);
+
+                indices[idx++] = (short)(FirstVertex+(DMA_GridSize+1)+1);
+                indices[idx++] = (short)(FirstVertex+(DMA_GridSize+1));
+                indices[idx++] = (short)(FirstVertex);
+            }
+            else
+            {
+                indices[idx++] = (short)(FirstVertex);
+                indices[idx++] = (short)(FirstVertex+1);
+                indices[idx++] = (short)(FirstVertex+(DMA_GridSize+1));
+
+                indices[idx++] = (short)(FirstVertex+1);
+                indices[idx++] = (short)(FirstVertex+(DMA_GridSize+1)+1);
+                indices[idx++] = (short)(FirstVertex+(DMA_GridSize+1));
+            }
+        }
+
+        return new DistortionMesh(vertices, vertices.length,
+                indices, indices.length);
+    }
+
+}
diff --git a/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/DistortionSpec.java b/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/DistortionSpec.java
new file mode 100644
index 0000000..f62a3c1
--- /dev/null
+++ b/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/DistortionSpec.java
@@ -0,0 +1,178 @@
+/**
+ * Copyright 2014 JogAmp Community. All rights reserved.
+ *
+ * Redistribution and use in source and binary forms, with or without modification, are
+ * permitted provided that the following conditions are met:
+ *
+ *    1. Redistributions of source code must retain the above copyright notice, this list of
+ *       conditions and the following disclaimer.
+ *
+ *    2. Redistributions in binary form must reproduce the above copyright notice, this list
+ *       of conditions and the following disclaimer in the documentation and/or other materials
+ *       provided with the distribution.
+ *
+ *    3. Compliance with Oculus VR RIFT SDK LICENSE (see below)
+ *
+ * THIS SOFTWARE IS PROVIDED BY JogAmp Community ``AS IS'' AND ANY EXPRESS OR IMPLIED
+ * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
+ * FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL JogAmp Community OR
+ * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
+ * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+ * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
+ * ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
+ * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ *
+ * The views and conclusions contained in the software and documentation are those of the
+ * authors and should not be interpreted as representing official policies, either expressed
+ * or implied, of JogAmp Community.
+ *
+ * XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
+ *
+ * This file contains mathematical equations, comments and algorithms
+ * used in the Oculus VR RIFT SDK 0.3.2.
+ *
+ * Due to unknown legal status, the 'requested' Copyright tag and disclaimer
+ * below has been added.
+ *
+ * XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
+ *
+ * Copyright © 2014 Oculus VR, Inc. All rights reserved
+ *
+ * Oculus VR, Inc. Software Development Kit License Agreement
+ *
+ * Human-Readable Summary*:
+ *
+ * You are Free to:
+ *
+ * Use, modify, and distribute the Oculus VR Rift SDK in source and binary
+ * form with your applications/software.
+ *
+ * With the Following Restrictions:
+ *
+ * You can only distribute or re-distribute the source code to LibOVR in
+ * whole, not in part.
+ *
+ * Modifications to the Oculus VR Rift SDK in source or binary form must
+ * be shared with Oculus VR.
+ *
+ * If your applications cause health and safety issues, you may lose your
+ * right to use the Oculus VR Rift SDK, including LibOVR.
+ *
+ * The Oculus VR Rift SDK may not be used to interface with unapproved commercial
+ * virtual reality mobile or non-mobile products or hardware.
+
+ * * - This human-readable Summary is not a license. It is simply a convenient
+ * reference for understanding the full Oculus VR Rift SDK License Agreement.
+ * The Summary is written as a user-friendly interface to the full Oculus VR Rift
+ * SDK License below. This Summary itself has no legal value, and its contents do
+ * not appear in the actual license.
+ *
+ * Full-length Legal Copy may be found at:
+ *   http://www.oculusvr.com/licenses/LICENSE-3.1
+ *   http://jogamp.org/git/?p=oculusvr-sdk.git;a=blob;f=LICENSE.txt;hb=HEAD
+ *   Or within this repository: oculusvr-sdk/LICENSE.txt
+ *
+ * THIS RIFT SDK AND ANY COMPONENT THEREOF IS PROVIDED BY OCULUS VR AND
+ * ITS CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING,
+ * BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
+ * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL OCULUS VR AS THE
+ * COPYRIGHT OWNER OR ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
+ * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+ * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,
+ * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+ * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS RIFT
+ * SDK OR THE RIFT SDK DERIVATIVES, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ */
+package jogamp.opengl.oculusvr.stereo.lense;
+
+import javax.media.nativewindow.util.DimensionImmutable;
+
+import jogamp.opengl.util.stereo.GenericStereoDevice;
+
+import com.jogamp.opengl.math.VectorUtil;
+
+public class DistortionSpec {
+    public DistortionSpec(final LensConfig lens) {
+        this.lens = lens;
+        this.pixelsPerTanAngleAtCenter = new float[2];
+        this.tanEyeAngleScale = new float[2];
+        this.lensCenter = new float[2];
+    }
+
+    private static final float[] quaterHalf = new float[] { 0.25f, 0.5f };
+
+    final LensConfig lens;
+
+    final float[] pixelsPerTanAngleAtCenter;
+    final float[] tanEyeAngleScale;
+    final float[] lensCenter;
+
+    public static DistortionSpec[] CalculateDistortionSpec (final GenericStereoDevice.Config deviceConfig, final float[] eyeReliefInMeters) {
+        final DistortionSpec[] result = new DistortionSpec[2];
+        /// FIXME: Add 'pluggable' lense configuration
+        final LensConfig[] lensConfig = LensConfig.GenerateLensConfigFromEyeRelief(eyeReliefInMeters, LensConfig.DistortionEquation.CatmullRom10);
+        result[0] = CalculateDistortionSpec (deviceConfig, 0, eyeReliefInMeters[0], lensConfig[0]);
+        result[1] = CalculateDistortionSpec (deviceConfig, 1, eyeReliefInMeters[1], lensConfig[1]);
+        return result;
+    }
+
+
+    private static DistortionSpec CalculateDistortionSpec (final GenericStereoDevice.Config deviceConfig, final int eyeName,
+                                                           final float eyeReliefInMeters, final LensConfig lensConfig) {
+        // From eye relief, IPD and device characteristics, we get the distortion mapping.
+        // This distortion does the following things:
+        // 1. It undoes the distortion that happens at the edges of the lens.
+        // 2. It maps the undistorted field into "retina" space.
+        // So the input is a pixel coordinate - the physical pixel on the display itself.
+        // The output is the real-world direction of the ray from this pixel as it comes out of the lens and hits the eye.
+        // However we typically think of rays "coming from" the eye, so the direction (TanAngleX,TanAngleY,1) is the direction
+        //      that the pixel appears to be in real-world space, where AngleX and AngleY are relative to the straight-ahead vector.
+        // If your renderer is a raytracer, you can use this vector directly (normalize as appropriate).
+        // However in standard rasterisers, we have rendered a 2D image and are putting it in front of the eye,
+        //      so we then need a mapping from this space to the [-1,1] UV coordinate space, which depends on exactly
+        //      where "in space" the app wants to put that rendertarget.
+        //      Where in space, and how large this rendertarget is, is completely up to the app and/or user,
+        //      though of course we can provide some useful hints.
+
+        // TODO: Use IPD and eye relief to modify distortion (i.e. non-radial component)
+        // TODO: cope with lenses that don't produce collimated light.
+        //       This means that IPD relative to the lens separation changes the light vergence,
+        //       and so we actually need to change where the image is displayed.
+        final DistortionSpec localDistortion = new DistortionSpec(lensConfig);
+
+        final DimensionImmutable resolutionInPixels = deviceConfig.surfaceSizeInPixels;
+
+        final float[] pixelsPerMeter = new float[] { resolutionInPixels.getWidth() / deviceConfig.screenSizeInMeters[0],
+                                                     resolutionInPixels.getHeight() / deviceConfig.screenSizeInMeters[1] };
+
+        VectorUtil.scaleVec2(localDistortion.pixelsPerTanAngleAtCenter, pixelsPerMeter, localDistortion.lens.MetersPerTanAngleAtCenter);
+
+        // Same thing, scaled to [-1,1] for each eye, rather than pixels.
+
+        VectorUtil.scaleVec2(localDistortion.tanEyeAngleScale, quaterHalf,
+                VectorUtil.divVec2(new float[2], deviceConfig.screenSizeInMeters, localDistortion.lens.MetersPerTanAngleAtCenter) );
+
+
+        // <--------------left eye------------------><-ScreenGapSizeInMeters-><--------------right eye----------------->
+        // <------------------------------------------ScreenSizeInMeters.Width----------------------------------------->
+        //                            <------------interpupillaryDistanceInMeters-------------->
+        // <--centerFromLeftInMeters->
+        //                            ^
+        //                      Center of lens
+
+        // Find the lens centers in scale of [-1,+1] (NDC) in left eye.
+        final float visibleWidthOfOneEye = 0.5f * ( deviceConfig.screenSizeInMeters[0] );
+        final float centerFromLeftInMeters = ( deviceConfig.screenSizeInMeters[0] - deviceConfig.interpupillaryDistanceInMeters ) * 0.5f;
+        localDistortion.lensCenter[0] = (     centerFromLeftInMeters / visibleWidthOfOneEye          ) * 2.0f - 1.0f;
+        localDistortion.lensCenter[1] = ( deviceConfig.pupilCenterFromScreenTopInMeters  / deviceConfig.screenSizeInMeters[1] ) * 2.0f - 1.0f;
+        if ( 1 == eyeName ) {
+            localDistortion.lensCenter[0] = -localDistortion.lensCenter[0];
+        }
+
+        return localDistortion;
+    }
+
+
+}
\ No newline at end of file
diff --git a/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/LensConfig.java b/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/LensConfig.java
new file mode 100644
index 0000000..f90470b
--- /dev/null
+++ b/jogl/src/classes/jogamp/opengl/oculusvr/stereo/lense/LensConfig.java
@@ -0,0 +1,600 @@
+/**
+ * Copyright 2014 JogAmp Community. All rights reserved.
+ *
+ * Redistribution and use in source and binary forms, with or without modification, are
+ * permitted provided that the following conditions are met:
+ *
+ *    1. Redistributions of source code must retain the above copyright notice, this list of
+ *       conditions and the following disclaimer.
+ *
+ *    2. Redistributions in binary form must reproduce the above copyright notice, this list
+ *       of conditions and the following disclaimer in the documentation and/or other materials
+ *       provided with the distribution.
+ *
+ *    3. Compliance with Oculus VR RIFT SDK LICENSE (see below)
+ *
+ * THIS SOFTWARE IS PROVIDED BY JogAmp Community ``AS IS'' AND ANY EXPRESS OR IMPLIED
+ * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
+ * FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL JogAmp Community OR
+ * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
+ * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+ * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
+ * ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
+ * ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ *
+ * The views and conclusions contained in the software and documentation are those of the
+ * authors and should not be interpreted as representing official policies, either expressed
+ * or implied, of JogAmp Community.
+ *
+ * XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
+ *
+ * This file contains mathematical equations, comments and algorithms
+ * used in the Oculus VR RIFT SDK 0.3.2.
+ *
+ * Due to unknown legal status, the 'requested' Copyright tag and disclaimer
+ * below has been added.
+ *
+ * XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
+ *
+ * Copyright © 2014 Oculus VR, Inc. All rights reserved
+ *
+ * Oculus VR, Inc. Software Development Kit License Agreement
+ *
+ * Human-Readable Summary*:
+ *
+ * You are Free to:
+ *
+ * Use, modify, and distribute the Oculus VR Rift SDK in source and binary
+ * form with your applications/software.
+ *
+ * With the Following Restrictions:
+ *
+ * You can only distribute or re-distribute the source code to LibOVR in
+ * whole, not in part.
+ *
+ * Modifications to the Oculus VR Rift SDK in source or binary form must
+ * be shared with Oculus VR.
+ *
+ * If your applications cause health and safety issues, you may lose your
+ * right to use the Oculus VR Rift SDK, including LibOVR.
+ *
+ * The Oculus VR Rift SDK may not be used to interface with unapproved commercial
+ * virtual reality mobile or non-mobile products or hardware.
+
+ * * - This human-readable Summary is not a license. It is simply a convenient
+ * reference for understanding the full Oculus VR Rift SDK License Agreement.
+ * The Summary is written as a user-friendly interface to the full Oculus VR Rift
+ * SDK License below. This Summary itself has no legal value, and its contents do
+ * not appear in the actual license.
+ *
+ * Full-length Legal Copy may be found at:
+ *   http://www.oculusvr.com/licenses/LICENSE-3.1
+ *   http://jogamp.org/git/?p=oculusvr-sdk.git;a=blob;f=LICENSE.txt;hb=HEAD
+ *   Or within this repository: oculusvr-sdk/LICENSE.txt
+ *
+ * THIS RIFT SDK AND ANY COMPONENT THEREOF IS PROVIDED BY OCULUS VR AND
+ * ITS CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING,
+ * BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
+ * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL OCULUS VR AS THE
+ * COPYRIGHT OWNER OR ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
+ * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+ * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,
+ * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+ * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+ * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS RIFT
+ * SDK OR THE RIFT SDK DERIVATIVES, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ */
+package jogamp.opengl.oculusvr.stereo.lense;
+
+public class LensConfig {
+    public static enum DistortionEquation { RecipPoly4, CatmullRom10 };
+
+    public static final int NumCoefficients = 11;
+
+    private DistortionEquation eqn;
+    /* pp */ float MetersPerTanAngleAtCenter;
+    private final float[] K;
+    private final float[] InvK;
+
+    private float MaxR;       // The highest R you're going to query for - the curve is unpredictable beyond it.
+    private float MaxInvR;
+
+    // Additional per-channel scaling is applied after distortion:
+    //  Index [0] - Red channel constant coefficient.
+    //  Index [1] - Red channel r^2 coefficient.
+    //  Index [2] - Blue channel constant coefficient.
+    //  Index [3] - Blue channel r^2 coefficient.
+    private final float[] ChromaticAberration = new float[4];
+
+    public LensConfig() {
+        this.K = new float[NumCoefficients];
+        this.InvK = new float[NumCoefficients];
+        SetToIdentity();
+    }
+
+    public LensConfig(final DistortionEquation eqn, final float MetersPerTanAngleAtCenter, final float[] K) {
+        this.K = new float[NumCoefficients];
+        this.InvK = new float[NumCoefficients];
+        SetToIdentity();
+        this.eqn = eqn;
+        this.MetersPerTanAngleAtCenter = MetersPerTanAngleAtCenter;
+        System.arraycopy(K, 0, this.K, 0, K.length);
+
+        // Chromatic aberration doesn't seem to change with eye relief.
+        ChromaticAberration[0]        = -0.006f;
+        ChromaticAberration[1]        =  0.0f;
+        ChromaticAberration[2]        =  0.014f;
+        ChromaticAberration[3]        =  0.0f;
+    }
+
+    private void SetUpInverseApprox() {
+        switch ( eqn )
+        {
+            case RecipPoly4: {
+                final float[] sampleR = new float[4];
+                final float[] sampleRSq = new float[4];
+                final float[] sampleInv = new float[4];
+                final float[] sampleFit = new float[4];
+                final float maxR = MaxInvR;
+
+                // Found heuristically...
+                sampleR[0] = 0.0f;
+                sampleR[1] = maxR * 0.4f;
+                sampleR[2] = maxR * 0.8f;
+                sampleR[3] = maxR * 1.5f;
+                for ( int i = 0; i < 4; i++ ) {
+                    sampleRSq[i] = sampleR[i] * sampleR[i];
+                    sampleInv[i] = DistortionFnInverse ( sampleR[i] );
+                    sampleFit[i] = sampleR[i] / sampleInv[i];
+                }
+                sampleFit[0] = 1.0f;
+                FitCubicPolynomial ( InvK, sampleRSq, sampleFit );
+            }
+            break;
+            case CatmullRom10: {
+                final int NumSegments = NumCoefficients;
+                for ( int i = 1; i < NumSegments; i++ ) {
+                    final float scaledRsq = i;
+                    final float rsq = scaledRsq * MaxInvR * MaxInvR / ( NumSegments - 1);
+                    final float r = (float)Math.sqrt ( rsq );
+                    final float inv = DistortionFnInverse ( r );
+                    InvK[i] = inv / r;
+                    InvK[0] = 1.0f;     // TODO: fix this.
+                }
+            }
+            break;
+        default:
+            throw new InternalError("unsupported EQ "+eqn);
+        }
+
+    }
+
+    private void SetToIdentity() {
+        for ( int i = 0; i < NumCoefficients; i++ )
+        {
+            K[i] = 0.0f;
+            InvK[i] = 0.0f;
+        }
+        eqn  = DistortionEquation.RecipPoly4;
+        K[0] = 1.0f;
+        InvK[0] = 1.0f;
+        MaxR = 1.0f;
+        MaxInvR = 1.0f;
+        ChromaticAberration[0] = 0.0f;
+        ChromaticAberration[1] = 0.0f;
+        ChromaticAberration[2] = 0.0f;
+        ChromaticAberration[3] = 0.0f;
+        MetersPerTanAngleAtCenter = 0.05f;
+    }
+
+    // DistortionFn applies distortion to the argument.
+    // Input: the distance in TanAngle/NIC space from the optical center to the input pixel.
+    // Output: the resulting distance after distortion.
+    public float DistortionFn(final float r) {
+        return r * DistortionFnScaleRadiusSquared ( r * r );
+    }
+
+    // The result is a scaling applied to the distance.
+    public float DistortionFnInverseApprox(final float r) {
+        final float rsq = r * r;
+        final float scale;
+
+        switch ( eqn )
+        {
+            case RecipPoly4: {
+                scale = 1.0f / ( InvK[0] + rsq * ( InvK[1] + rsq * ( InvK[2] + rsq * InvK[3] ) ) );
+            }
+            break;
+            case CatmullRom10: {
+                // A Catmull-Rom spline through the values 1.0, K[1], K[2] ... K[9]
+                // evenly spaced in R^2 from 0.0 to MaxR^2
+                // K[0] controls the slope at radius=0.0, rather than the actual value.
+                final int NumSegments = NumCoefficients;
+                final float scaledRsq = (NumSegments-1) * rsq / ( MaxInvR * MaxInvR );
+                scale = EvalCatmullRom10Spline ( InvK, scaledRsq );
+            }
+            break;
+        default:
+            throw new InternalError("unsupported EQ "+eqn);
+        }
+        return r * scale;
+    }
+
+    // DistortionFnInverse computes the inverse of the distortion function on an argument.
+    public float DistortionFnInverse(final float r) {
+        float s, d;
+        float delta = r * 0.25f;
+
+        // Better to start guessing too low & take longer to converge than too high
+        // and hit singularities. Empirically, r * 0.5f is too high in some cases.
+        s = r * 0.25f;
+        d = Math.abs(r - DistortionFn(s));
+
+        for (int i = 0; i < 20; i++)
+        {
+            final float sUp   = s + delta;
+            final float sDown = s - delta;
+            final float dUp   = Math.abs(r - DistortionFn(sUp));
+            final float dDown = Math.abs(r - DistortionFn(sDown));
+
+            if (dUp < d)
+            {
+                s = sUp;
+                d = dUp;
+            }
+            else if (dDown < d)
+            {
+                s = sDown;
+                d = dDown;
+            }
+            else
+            {
+                delta *= 0.5f;
+            }
+        }
+        return s;
+    }
+
+    // The result is a scaling applied to the distance from the center of the lens.
+    public float DistortionFnScaleRadiusSquared (final float rsq) {
+        final float scale;
+        switch ( eqn )
+        {
+            case RecipPoly4: {
+                scale = 1.0f / ( K[0] + rsq * ( K[1] + rsq * ( K[2] + rsq * K[3] ) ) );
+            }
+            break;
+            case CatmullRom10: {
+                // A Catmull-Rom spline through the values 1.0, K[1], K[2] ... K[10]
+                // evenly spaced in R^2 from 0.0 to MaxR^2
+                // K[0] controls the slope at radius=0.0, rather than the actual value.
+                final int NumSegments = LensConfig.NumCoefficients;
+                final float scaledRsq = (NumSegments-1) * rsq / ( MaxR * MaxR );
+                scale = EvalCatmullRom10Spline ( K, scaledRsq );
+            }
+            break;
+        default:
+            throw new InternalError("unsupported EQ "+eqn);
+        }
+        return scale;
+    }
+
+    // x,y,z components map to r,g,b
+    public float[] DistortionFnScaleRadiusSquaredChroma(final float radiusSquared) {
+        final float scale = DistortionFnScaleRadiusSquared ( radiusSquared );
+        final float[] scaleRGB = new float[3];
+        scaleRGB[0] = scale * ( 1.0f + ChromaticAberration[0] + radiusSquared * ChromaticAberration[1] );     // Red
+        scaleRGB[1] = scale;                                                                                  // Green
+        scaleRGB[2] = scale * ( 1.0f + ChromaticAberration[2] + radiusSquared * ChromaticAberration[3] );     // Blue
+        return scaleRGB;
+    }
+
+    /**
+    // Inputs are 4 points (pFitX[0],pFitY[0]) through (pFitX[3],pFitY[3])
+    // Result is four coefficients in pResults[0] through pResults[3] such that
+    //      y = pResult[0] + x * ( pResult[1] + x * ( pResult[2] + x * ( pResult[3] ) ) );
+    // passes through all four input points.
+    // Return is true if it succeeded, false if it failed (because two control points
+    // have the same pFitX value).
+     *
+     * @param pResult
+     * @param pFitX
+     * @param pFitY
+     * @return
+     */
+    private static boolean FitCubicPolynomial ( final float[/*4*/] pResult, final float[/*4*/] pFitX, final float[/*4*/] pFitY ) {
+        final float d0 = ( ( pFitX[0]-pFitX[1] ) * ( pFitX[0]-pFitX[2] ) * ( pFitX[0]-pFitX[3] ) );
+        final float d1 = ( ( pFitX[1]-pFitX[2] ) * ( pFitX[1]-pFitX[3] ) * ( pFitX[1]-pFitX[0] ) );
+        final float d2 = ( ( pFitX[2]-pFitX[3] ) * ( pFitX[2]-pFitX[0] ) * ( pFitX[2]-pFitX[1] ) );
+        final float d3 = ( ( pFitX[3]-pFitX[0] ) * ( pFitX[3]-pFitX[1] ) * ( pFitX[3]-pFitX[2] ) );
+
+        if ( ( d0 == 0.0f ) || ( d1 == 0.0f ) || ( d2 == 0.0f ) || ( d3 == 0.0f ) )
+        {
+            return false;
+        }
+
+        final float f0 = pFitY[0] / d0;
+        final float f1 = pFitY[1] / d1;
+        final float f2 = pFitY[2] / d2;
+        final float f3 = pFitY[3] / d3;
+
+        pResult[0] = -( f0*pFitX[1]*pFitX[2]*pFitX[3]
+                      + f1*pFitX[0]*pFitX[2]*pFitX[3]
+                      + f2*pFitX[0]*pFitX[1]*pFitX[3]
+                      + f3*pFitX[0]*pFitX[1]*pFitX[2] );
+        pResult[1] = f0*(pFitX[1]*pFitX[2] + pFitX[2]*pFitX[3] + pFitX[3]*pFitX[1])
+                   + f1*(pFitX[0]*pFitX[2] + pFitX[2]*pFitX[3] + pFitX[3]*pFitX[0])
+                   + f2*(pFitX[0]*pFitX[1] + pFitX[1]*pFitX[3] + pFitX[3]*pFitX[0])
+                   + f3*(pFitX[0]*pFitX[1] + pFitX[1]*pFitX[2] + pFitX[2]*pFitX[0]);
+        pResult[2] = -( f0*(pFitX[1]+pFitX[2]+pFitX[3])
+                      + f1*(pFitX[0]+pFitX[2]+pFitX[3])
+                      + f2*(pFitX[0]+pFitX[1]+pFitX[3])
+                      + f3*(pFitX[0]+pFitX[1]+pFitX[2]) );
+        pResult[3] = f0 + f1 + f2 + f3;
+
+        return true;
+    }
+
+    private static float EvalCatmullRom10Spline ( final float[] K, final float scaledVal ) {
+        final int NumSegments = NumCoefficients;
+
+        float scaledValFloor = (float)Math.floor( scaledVal );
+        scaledValFloor = Math.max( 0.0f, Math.min( NumSegments-1, scaledValFloor ) );
+        final float t = scaledVal - scaledValFloor;
+        final int k = (int)scaledValFloor;
+
+        float p0, p1;
+        float m0, m1;
+        switch ( k )
+        {
+        case 0:
+            // Curve starts at 1.0 with gradient K[1]-K[0]
+            p0 = 1.0f;
+            m0 =        ( K[1] - K[0] );    // general case would have been (K[1]-K[-1])/2
+            p1 = K[1];
+            m1 = 0.5f * ( K[2] - K[0] );
+            break;
+        default:
+            // General case
+            p0 = K[k  ];
+            m0 = 0.5f * ( K[k+1] - K[k-1] );
+            p1 = K[k+1];
+            m1 = 0.5f * ( K[k+2] - K[k  ] );
+            break;
+        case NumSegments-2:
+            // Last tangent is just the slope of the last two points.
+            p0 = K[NumSegments-2];
+            m0 = 0.5f * ( K[NumSegments-1] - K[NumSegments-2] );
+            p1 = K[NumSegments-1];
+            m1 = K[NumSegments-1] - K[NumSegments-2];
+            break;
+        case NumSegments-1:
+            // Beyond the last segment it's just a straight line
+            p0 = K[NumSegments-1];
+            m0 = K[NumSegments-1] - K[NumSegments-2];
+            p1 = p0 + m0;
+            m1 = m0;
+            break;
+        }
+
+        final float omt = 1.0f - t;
+        final float res  = ( p0 * ( 1.0f + 2.0f *   t ) + m0 *   t ) * omt * omt
+                + ( p1 * ( 1.0f + 2.0f * omt ) - m1 * omt ) *   t *   t;
+
+        return res;
+    }
+
+    /** FIXME: Add 'pluggable' lense configuration */
+    public static LensConfig[] GenerateLensConfigFromEyeRelief(final float[] eyeReliefInMeters, final DistortionEquation eqn) {
+        final LensConfig[] result = new LensConfig[2];
+        final DistortionDescriptor[] distortions = LensConfig.CreateDistortionDescriptorsforOVRDK1_CupsABC();
+        result[0] = GenerateLensConfigFromEyeRelief(eyeReliefInMeters[0], distortions, eqn);
+        result[1] = GenerateLensConfigFromEyeRelief(eyeReliefInMeters[1], distortions, eqn);
+        return result;
+    }
+
+    private static LensConfig GenerateLensConfigFromEyeRelief(final float eyeReliefInMeters, final DistortionDescriptor[] distortions, final DistortionEquation eqn) {
+        final int numDistortions = distortions.length;
+        final int defaultDistortion = 0; // index of the default distortion curve to use if zero eye relief supplied
+
+        DistortionDescriptor pUpper = null;
+        DistortionDescriptor pLower = null;
+        float lerpVal = 0.0f;
+        if (eyeReliefInMeters == 0)
+        {   // Use a constant default distortion if an invalid eye-relief is supplied
+            pLower = distortions[defaultDistortion];
+            pUpper = distortions[defaultDistortion];
+            lerpVal = 0.0f;
+        } else {
+            for ( int i = 0; i < numDistortions-1; i++ )
+            {
+                assert( distortions[i].eyeRelief < distortions[i+1].eyeRelief );
+                if ( ( distortions[i].eyeRelief <= eyeReliefInMeters ) && ( distortions[i+1].eyeRelief > eyeReliefInMeters ) )
+                {
+                    pLower = distortions[i];
+                    pUpper = distortions[i+1];
+                    lerpVal = ( eyeReliefInMeters - pLower.eyeRelief ) / ( pUpper.eyeRelief - pLower.eyeRelief );
+                    // No break here - I want the ASSERT to check everything every time!
+                }
+            }
+        }
+
+        if ( pUpper == null )
+        {
+            // Do not extrapolate, just clamp - slightly worried about people putting in bogus settings.
+            if ( distortions[0].eyeRelief > eyeReliefInMeters )
+            {
+                pLower = distortions[0];
+                pUpper = distortions[0];
+            }
+            else
+            {
+                assert ( distortions[numDistortions-1].eyeRelief <= eyeReliefInMeters );
+                pLower = distortions[numDistortions-1];
+                pUpper = distortions[numDistortions-1];
+            }
+            lerpVal = 0.0f;
+        }
+        final float invLerpVal = 1.0f - lerpVal;
+
+        pLower.config.MaxR = pLower.maxRadius;
+        pUpper.config.MaxR = pUpper.maxRadius;
+
+        final LensConfig result = new LensConfig();
+        // Where is the edge of the lens - no point modelling further than this.
+        final float maxValidRadius = invLerpVal * pLower.maxRadius + lerpVal * pUpper.maxRadius;
+        result.MaxR = maxValidRadius;
+
+        switch ( eqn )
+        {
+            case RecipPoly4:{
+                    // Lerp control points and fit an equation to them.
+                    final float[] fitX = new float[4];
+                    final float[] fitY = new float[4];
+                    fitX[0] = 0.0f;
+                    fitY[0] = 1.0f;
+                    for ( int ctrlPt = 1; ctrlPt < 4; ctrlPt ++ )
+                    {
+                        final float radiusLerp = invLerpVal * pLower.sampleRadius[ctrlPt-1] + lerpVal * pUpper.sampleRadius[ctrlPt-1];
+                        final float radiusLerpSq = radiusLerp * radiusLerp;
+                        final float fitYLower = pLower.config.DistortionFnScaleRadiusSquared ( radiusLerpSq );
+                        final float fitYUpper = pUpper.config.DistortionFnScaleRadiusSquared ( radiusLerpSq );
+                        fitX[ctrlPt] = radiusLerpSq;
+                        fitY[ctrlPt] = 1.0f / ( invLerpVal * fitYLower + lerpVal * fitYUpper );
+                    }
+                    result.eqn = DistortionEquation.RecipPoly4;
+                    final boolean bSuccess = LensConfig.FitCubicPolynomial ( result.K, fitX, fitY );
+                    assert ( bSuccess );
+
+                    // Set up the fast inverse.
+                    final float maxRDist = result.DistortionFn ( maxValidRadius );
+                    result.MaxInvR = maxRDist;
+                    result.SetUpInverseApprox();
+                }
+                break;
+            case CatmullRom10: {
+                // Evenly sample & lerp points on the curve.
+                final int NumSegments = LensConfig.NumCoefficients;
+                result.MaxR = maxValidRadius;
+                // Directly interpolate the K0 values
+                result.K[0] = invLerpVal * pLower.config.K[0] + lerpVal * pUpper.config.K[0];
+
+                // Sample and interpolate the distortion curves to derive K[1] ... K[n]
+                for ( int ctrlPt = 1; ctrlPt < NumSegments; ctrlPt++ )
+                {
+                    final float radiusSq = ( (float)ctrlPt / (float)(NumSegments-1) ) * maxValidRadius * maxValidRadius;
+                    final float fitYLower = pLower.config.DistortionFnScaleRadiusSquared ( radiusSq );
+                    final float fitYUpper = pUpper.config.DistortionFnScaleRadiusSquared ( radiusSq );
+                    final float fitLerp = invLerpVal * fitYLower + lerpVal * fitYUpper;
+                    result.K[ctrlPt] = fitLerp;
+                }
+
+                result.eqn = DistortionEquation.CatmullRom10;
+
+                for ( int ctrlPt = 1; ctrlPt < NumSegments; ctrlPt++ )
+                {
+                    final float radiusSq = ( (float)ctrlPt / (float)(NumSegments-1) ) * maxValidRadius * maxValidRadius;
+                    final float val = result.DistortionFnScaleRadiusSquared ( radiusSq );
+                    assert ( Math.abs( val - result.K[ctrlPt] ) < 0.0001f );
+                }
+
+                // Set up the fast inverse.
+                final float maxRDist = result.DistortionFn ( maxValidRadius );
+                result.MaxInvR = maxRDist;
+                result.SetUpInverseApprox();
+            }
+            break;
+        default:
+            throw new InternalError("unsupported EQ "+eqn);
+        }
+
+        // Chromatic aberration.
+        result.ChromaticAberration[0] = invLerpVal * pLower.config.ChromaticAberration[0] + lerpVal * pUpper.config.ChromaticAberration[0];
+        result.ChromaticAberration[1] = invLerpVal * pLower.config.ChromaticAberration[1] + lerpVal * pUpper.config.ChromaticAberration[1];
+        result.ChromaticAberration[2] = invLerpVal * pLower.config.ChromaticAberration[2] + lerpVal * pUpper.config.ChromaticAberration[2];
+        result.ChromaticAberration[3] = invLerpVal * pLower.config.ChromaticAberration[3] + lerpVal * pUpper.config.ChromaticAberration[3];
+
+        // Scale.
+        result.MetersPerTanAngleAtCenter =  pLower.config.MetersPerTanAngleAtCenter * invLerpVal +
+                                            pUpper.config.MetersPerTanAngleAtCenter * lerpVal;
+
+        return result;
+    }
+
+    public static class DistortionDescriptor {
+        public DistortionDescriptor(final LensConfig lens, final float eyeRelief,
+                                    final float[] sampleRadius, final float maxRadius) {
+            this.config = lens;
+            this.eyeRelief = eyeRelief;
+            this.sampleRadius = sampleRadius;
+            this.maxRadius = maxRadius;
+        }
+
+        final LensConfig config;
+        final float eyeRelief;
+        final float[] sampleRadius;
+        final float maxRadius;
+    }
+
+    /*** Hardcoded OculusVR DK1 A, B, C eye cups (lenses) */
+    public static DistortionDescriptor[] CreateDistortionDescriptorsforOVRDK1_CupsABC() {
+        return new DistortionDescriptor[] {
+            // Tuned at minimum dial setting - extended to r^2 == 1.8
+            new DistortionDescriptor(
+                    new LensConfig(DistortionEquation.CatmullRom10,
+                                   0.0425f, // MetersPerTanAngleAtCenter
+                                   new float[] { 1.0000f,   // K00
+                                                 1.06505f,  // K01
+                                                 1.14725f,  // K02
+                                                 1.2705f,   // K03
+                                                 1.48f,     // K04
+                                                 1.87f,     // K05
+                                                 2.534f,    // K06
+                                                 3.6f,      // K07
+                                                 5.1f,      // K08
+                                                 7.4f,      // K09
+                                                 11.0f} ),  // K10
+                    0.012760465f - 0.005f, // eyeRelief
+                    new float[] { 0.222717149f, 0.512249443f, 0.712694878f }, // sampleRadius
+                    (float)Math.sqrt(1.8f) ), // maxRadius
+            // Tuned at middle dial setting
+            new DistortionDescriptor(
+                    new LensConfig(DistortionEquation.CatmullRom10,
+                                   0.0425f, // MetersPerTanAngleAtCenter
+                                   new float[] { 1.0000f,       // K00
+                                                 1.032407264f,  // K01
+                                                 1.07160462f,   // K02
+                                                 1.11998388f,   // K03
+                                                 1.1808606f,    // K04
+                                                 1.2590494f,    // K05
+                                                 1.361915f,     // K06
+                                                 1.5014339f,    // K07
+                                                 1.6986004f,    // K08
+                                                 1.9940577f,    // K09
+                                                 2.4783147f} ), // K10
+                    0.012760465f, // eyeRelief
+                    new float[] { 0.222717149f, 0.512249443f, 0.712694878f }, // sampleRadius
+                    1.0f ), // maxRadius
+            // Tuned at maximum dial setting
+            new DistortionDescriptor(
+                    new LensConfig(DistortionEquation.CatmullRom10,
+                                   0.0425f, // MetersPerTanAngleAtCenter
+                                   new float[] { 1.0102f,       // K00
+                                                 1.0371f,  // K01
+                                                 1.0831f,   // K02
+                                                 1.1353f,   // K03
+                                                 1.2f,    // K04
+                                                 1.2851f,    // K05
+                                                 1.3979f,     // K06
+                                                 1.56f,    // K07
+                                                 1.8f,    // K08
+                                                 2.25f,    // K09
+                                                 3.0f} ), // K10
+                    0.012760465f + 0.005f, // eyeRelief
+                    new float[] { 0.222717149f, 0.512249443f, 0.712694878f }, // sampleRadius
+                    1.0f ), // maxRadius
+        };
+    }
+
+}
\ No newline at end of file