Fix inverse Lut1D optimization bug (#1743)

Signed-off-by: Doug Walker <doug.walker@autodesk.com>

Signed-off-by: Doug Walker <doug.walker@autodesk.com>

src/OpenColorIO/OpOptimizers.cpp

@@ -12,6 +12,7 @@
 #include "Op.h"
 #include "ops/lut1d/Lut1DOp.h"
 #include "ops/lut3d/Lut3DOp.h"
+#include "ops/range/RangeOp.h"
 
 namespace OCIO_NAMESPACE
 {
@@ -241,7 +242,38 @@ int RemoveInverseOps(OpRcPtrVec & opVec, OptimizationFlags oFlags)
             // When a pair of inverse ops is removed, we want the optimized ops to give the
             // same result as the original.  For certain ops such as Lut1D or Log this may
             // mean inserting a Range to emulate the clamping done by the original ops.
-            auto replacedBy = op1->getIdentityReplacement();
+
+            OpRcPtr replacedBy;
+            if (type1 == OpData::Lut1DType)
+            {
+                // Lut1D gets special handling so that both halfs of the pair are available.
+                // Only the inverse LUT has the values needed to generate the replacement.
+
+                ConstLut1DOpDataRcPtr lut1 = OCIO_DYNAMIC_POINTER_CAST<const Lut1DOpData>(op1->data());
+                ConstLut1DOpDataRcPtr lut2 = OCIO_DYNAMIC_POINTER_CAST<const Lut1DOpData>(op2->data());
+
+                OpDataRcPtr opData = lut1->getPairIdentityReplacement(lut2);
+
+                OpRcPtrVec ops;
+                if (opData->getType() == OpData::MatrixType)
+                {
+                    // No-op that will be optimized.
+                    auto mat = OCIO_DYNAMIC_POINTER_CAST<MatrixOpData>(opData);
+                    CreateMatrixOp(ops, mat, TRANSFORM_DIR_FORWARD);
+                }
+                else if (opData->getType() == OpData::RangeType)
+                {
+                    // Clamping op.
+                    auto range = OCIO_DYNAMIC_POINTER_CAST<RangeOpData>(opData);
+                    CreateRangeOp(ops, range, TRANSFORM_DIR_FORWARD);
+                }
+                replacedBy = ops[0];
+            }
+            else
+            {
+                replacedBy = op1->getIdentityReplacement();
+            }
+
             replacedBy->finalize();
             if (replacedBy->isNoOp())
             {

src/OpenColorIO/ops/lut1d/Lut1DOpData.cpp

@@ -279,6 +279,86 @@ OpDataRcPtr Lut1DOpData::getIdentityReplacement() const
     return res;
 }
 
+OpDataRcPtr Lut1DOpData::getPairIdentityReplacement(ConstLut1DOpDataRcPtr & lut2) const
+{
+    OpDataRcPtr res;
+    if (isInputHalfDomain())
+    {
+        // TODO: If a half-domain LUT has a flat spot, it would be more appropriate
+        // to use a Range, since some areas would be clamped in a round-trip.
+        // Currently leaving this a Matrix since it is a potential work-around
+        // for situations where you want a pair identity of LUTs to be totally
+        // removed, even if it omits some clamping at extreme values.
+        res = std::make_shared<MatrixOpData>();
+    }
+    else
+    {
+        // Note that the ops have been finalized by the time this is called,
+        // Therefore, for an inverse Lut1D, it means initializeFromForward() has
+        // been called and so any reversals have been converted to flat regions.
+        // Therefore, the first and last LUT entries are the extreme values and
+        // the ComponentProperties has been initialized, but only for the op
+        // whose direction is INVERSE.
+        const Lut1DOpData * invLut = (m_direction == TRANSFORM_DIR_INVERSE)
+            ? this: lut2.get();
+        const ComponentProperties & redProperties = invLut->getRedProperties();
+        const unsigned long length = invLut->getArray().getLength();
+
+        // If the start or end of the LUT contains a flat region, that will cause
+        // a round-trip to be limited.
+
+        double minValue = 0.;
+        double maxValue = 1.;
+        switch (m_direction)
+        {
+        case TRANSFORM_DIR_FORWARD: // Fwd Lut1D -> Inv Lut1D
+        {
+            // A round-trip in this order will impose at least a clamp to [0,1]
+            // based on what happens entering the first Fwd Lut1D.  However, the
+            // clamping may be to an even narrower range if there are flat regions.
+            //
+            // The flat region limitation is imposed based on the where it falls
+            // relative to the [0,1] input domain.
+
+            // TODO: A RangeOp has one min & max for all channels, whereas a Lut1D may 
+            // have three independent channels.  Potentially could look at all chans
+            // and take the extrema of each?  For now, just using the first channel.
+            const unsigned long minIndex = redProperties.startDomain;
+            const unsigned long maxIndex = redProperties.endDomain;
+
+            minValue = (double)minIndex / (length - 1);
+            maxValue = (double)maxIndex / (length - 1);
+            break;
+        }
+        case TRANSFORM_DIR_INVERSE: // Inv Lut1D -> Fwd Lut1D
+        {
+            // A round-trip in this order will impose a clamp, but it may be to
+            // bounds outside of [0,1] since the Fwd LUT may contain values outside
+            // [0,1] and so the Inv LUT will accept inputs on that extended range.
+            //
+            // The flat region limitation is imposed based on the output range.
+
+            const bool isIncreasing = redProperties.isIncreasing;
+            const unsigned long maxChannels = invLut->getArray().getMaxColorComponents();
+            const unsigned long lastValIndex = (length - 1) * maxChannels;
+            // Note that the array for the invLut has had initializeFromForward()
+            // done and so any reversals have been converted to flat regions and
+            // the extrema are at the beginning & end of the LUT.
+            const Array::Values & lutValues = invLut->getArray().getValues();
+
+            // TODO: Currently only basing this on the red channel.
+            minValue = isIncreasing ? lutValues[0] : lutValues[lastValIndex];
+            maxValue = isIncreasing ? lutValues[lastValIndex] : lutValues[0];
+            break;
+        }
+        }
+
+        res = std::make_shared<RangeOpData>(minValue, maxValue,
+                                            minValue, maxValue);
+    }
+    return res;
+}
+
 void Lut1DOpData::setInputHalfDomain(bool isHalfDomain) noexcept
 {
     m_halfFlags = (isHalfDomain) ?

src/OpenColorIO/ops/lut1d/Lut1DOpData.h

@@ -179,6 +179,8 @@ class Lut1DOpData : public OpData
 
     OpDataRcPtr getIdentityReplacement() const override;
 
+    OpDataRcPtr getPairIdentityReplacement(ConstLut1DOpDataRcPtr & lut2) const;
+
     inline const ComponentProperties & getRedProperties() const
     {
         return m_componentProperties[0];

tests/cpu/OpOptimizers_tests.cpp

@@ -645,6 +645,63 @@ OCIO_ADD_TEST(OpOptimizers, lut1d_identity_replacement)
     }
 }
 
+OCIO_ADD_TEST(OpOptimizers, lut1d_identity_replacement_order)
+{
+    // See issue #1737, https://github.com/AcademySoftwareFoundation/OpenColorIO/issues/1737.
+
+    // This CTF contains a single LUT1D, inverse direction, normal (not half) domain.
+    // It contains values from -6 to +3.4.
+    const std::string fileName("lut1d_inverse_gpu.ctf");
+    OCIO::ContextRcPtr context = OCIO::Context::Create();
+
+    OCIO::OpRcPtrVec inv_ops;
+    OCIO_CHECK_NO_THROW(OCIO::BuildOpsTest(inv_ops, fileName, context,
+                                        // FWD direction simply means don't swap the direction, the
+                                        // file contains an inverse LUT1D and leave it that way.
+                                           OCIO::TRANSFORM_DIR_FORWARD));
+    OCIO::OpRcPtrVec fwd_ops;
+    OCIO_CHECK_NO_THROW(OCIO::BuildOpsTest(fwd_ops, fileName, context,
+                                           OCIO::TRANSFORM_DIR_INVERSE));
+
+    // Check forward LUT1D followed by inverse LUT1D.
+    {
+        OCIO::OpRcPtrVec fwd_inv_ops = fwd_ops;
+        fwd_inv_ops += inv_ops;
+
+        OCIO_CHECK_NO_THROW(fwd_inv_ops.finalize());
+        OCIO_CHECK_NO_THROW(fwd_inv_ops.optimize(OCIO::OPTIMIZATION_NONE));
+        OCIO_CHECK_EQUAL(fwd_inv_ops.size(), 2);  // no optmization was done
+
+        OCIO::OpRcPtrVec optOps = fwd_inv_ops.clone();
+        OCIO_CHECK_NO_THROW(optOps.finalize());
+        OCIO_CHECK_NO_THROW(optOps.optimize(OCIO::OPTIMIZATION_DEFAULT));
+        OCIO_CHECK_EQUAL(optOps.size(), 1);
+        OCIO_CHECK_EQUAL(optOps[0]->getInfo(), "<RangeOp>");
+
+        // Compare renders.
+        CompareRender(fwd_inv_ops, optOps, __LINE__, 1e-6f);
+    }
+
+    // Check inverse LUT1D followed by forward LUT1D.
+    {
+        OCIO::OpRcPtrVec inv_fwd_ops = inv_ops;
+        inv_fwd_ops += fwd_ops;
+
+        OCIO_CHECK_NO_THROW(inv_fwd_ops.finalize());
+        OCIO_CHECK_NO_THROW(inv_fwd_ops.optimize(OCIO::OPTIMIZATION_NONE));
+        OCIO_CHECK_EQUAL(inv_fwd_ops.size(), 2);  // no optmization was done
+
+        OCIO::OpRcPtrVec optOps = inv_fwd_ops.clone();
+        OCIO_CHECK_NO_THROW(optOps.finalize());
+        OCIO_CHECK_NO_THROW(optOps.optimize(OCIO::OPTIMIZATION_DEFAULT));
+        OCIO_CHECK_EQUAL(optOps.size(), 1);
+        OCIO_CHECK_EQUAL(optOps[0]->getInfo(), "<RangeOp>");
+
+        // Compare renders.
+        CompareRender(inv_fwd_ops, optOps, __LINE__, 1e-6f);
+    }
+}
+
 OCIO_ADD_TEST(OpOptimizers, lut1d_half_domain_keep_prior_range)
 {
     // A half-domain LUT should not allow removal of a prior range op.