fix(png): Improve png write with alpha is low (AcademySoftwareFoundation#3985)

We found that PNG output (which is required to be unassociated alpha)
had precision problems with low alpha values because it quantized to the
final bit depth before doing the alpha deassociation.

Explanation:

PNG requires associated alpha to be stored in the file. So if you have
(unassociated) pixel value [0.00235, 0.00106, 0.00117, 0.0025], that
means that we need to divide the colors by alpha before storing in the
file, so that's [0.94, 0.424, 0.468, 0.0025] if you have full precision.
But you don't, because PNG only supports uint8 and uint16. For uint8,
this should end up as [240, 108, 119, 1].

But it turns out that our png writer does the integer conversion before
the un-premultiplication, so actually we first quantize from [0.00235,
0.00106, 0.00117, 0.0025] to [1 0 0 1], and THEN we un-premultiply,
giving us the [255, 0, 0, 1] that made you think "Why is that red value
clipped? And where did green and blue go?"

The solution, then, is to deassociate as floats first, then quantize to
the final integer data type.

Also add a test for this, as well as a test that verifies correctly
writing of 16 bit files that I had some problems with in the last png
PR.

Signed-off-by: Larry Gritz <[email protected]>

Author: lgritz

src/png.imageio/pngoutput.cpp

@@ -89,8 +89,8 @@ class PNGOutput final : public ImageOutput {
     }
 
     template<class T>
-    void deassociateAlpha(T* data, int size, int channels, int alpha_channel,
-                          float gamma);
+    void deassociateAlpha(T* data, size_t npixels, int channels,
+                          int alpha_channel, float gamma);
 };
 
 
@@ -269,30 +269,33 @@ PNGOutput::close()
 
 template<class T>
 void
-PNGOutput::deassociateAlpha(T* data, int size, int channels, int alpha_channel,
-                            float gamma)
+PNGOutput::deassociateAlpha(T* data, size_t npixels, int channels,
+                            int alpha_channel, float gamma)
 {
-    unsigned int max = std::numeric_limits<T>::max();
     if (gamma == 1) {
-        for (int x = 0; x < size; ++x, data += channels)
-            if (data[alpha_channel])
-                for (int c = 0; c < channels; c++)
-                    if (c != alpha_channel) {
-                        unsigned int f = data[c];
-                        f              = (f * max) / data[alpha_channel];
-                        data[c]        = (T)std::min(max, f);
-                    }
+        for (size_t x = 0; x < npixels; ++x, data += channels) {
+            DataArrayProxy<T, float> val(data);
+            float alpha = val[alpha_channel];
+            if (alpha != 0.0f && alpha != 1.0f) {
+                for (int c = 0; c < channels; c++) {
+                    if (c != alpha_channel)
+                        val[c] = data[c] / alpha;
+                }
+            }
+        }
     } else {
-        for (int x = 0; x < size; ++x, data += channels)
-            if (data[alpha_channel]) {
+        for (size_t x = 0; x < npixels; ++x, data += channels) {
+            DataArrayProxy<T, float> val(data);
+            float alpha = val[alpha_channel];
+            if (alpha != 0.0f && alpha != 1.0f) {
                 // See associateAlpha() for an explanation.
-                float alpha_deassociate = pow((float)max / data[alpha_channel],
-                                              gamma);
-                for (int c = 0; c < channels; c++)
+                float alpha_deassociate = pow(1.0f / val[alpha_channel], gamma);
+                for (int c = 0; c < channels; c++) {
                     if (c != alpha_channel)
-                        data[c] = static_cast<T>(std::min(
-                            max, (unsigned int)(data[c] * alpha_deassociate)));
+                        val[c] = val[c] * alpha_deassociate;
+                }
             }
+        }
     }
 }
 
@@ -302,26 +305,42 @@ bool
 PNGOutput::write_scanline(int y, int z, TypeDesc format, const void* data,
                           stride_t xstride)
 {
-    y -= m_spec.y;
     m_spec.auto_stride(xstride, format, spec().nchannels);
     const void* origdata = data;
+    if (format == TypeUnknown)
+        format = m_spec.format;
+
+    // PNG specifically dictates unassociated (un-"premultiplied") alpha.
+    // If we need to unassociate alpha, do it in float.
+    std::unique_ptr<float[]> unassoc_scratch;
+    if (m_convert_alpha) {
+        size_t nvals     = size_t(m_spec.width) * size_t(m_spec.nchannels);
+        float* floatvals = nullptr;
+        if (nvals * sizeof(float) <= (1 << 16)) {
+            floatvals = OIIO_ALLOCA(float, nvals);  // small enough for stack
+        } else {
+            unassoc_scratch.reset(new float[nvals]);
+            floatvals = unassoc_scratch.get();
+        }
+        // Contiguize and convert to float
+        OIIO::convert_image(m_spec.nchannels, m_spec.width, 1, 1, data, format,
+                            xstride, AutoStride, AutoStride, floatvals,
+                            TypeFloat, AutoStride, AutoStride, AutoStride);
+        // Deassociate alpha
+        deassociateAlpha(floatvals, size_t(m_spec.width), m_spec.nchannels,
+                         m_spec.alpha_channel, m_gamma);
+        data    = floatvals;
+        format  = TypeFloat;
+        xstride = size_t(m_spec.nchannels) * sizeof(float);
+    }
+
     data = to_native_scanline(format, data, xstride, m_scratch, m_dither, y, z);
     if (data == origdata && (m_convert_alpha || m_need_swap)) {
         m_scratch.assign((unsigned char*)data,
                          (unsigned char*)data + m_spec.scanline_bytes());
         data = &m_scratch[0];
     }
 
-    // PNG specifically dictates unassociated (un-"premultiplied") alpha
-    if (m_convert_alpha) {
-        if (m_spec.format == TypeDesc::UINT16)
-            deassociateAlpha((unsigned short*)data, m_spec.width,
-                             m_spec.nchannels, m_spec.alpha_channel, m_gamma);
-        else
-            deassociateAlpha((unsigned char*)data, m_spec.width,
-                             m_spec.nchannels, m_spec.alpha_channel, m_gamma);
-    }
-
     // PNG is always big endian
     if (m_need_swap)
         swap_endian((unsigned short*)data, m_spec.width * m_spec.nchannels);
@@ -354,27 +373,41 @@ PNGOutput::write_scanlines(int ybegin, int yend, int z, TypeDesc format,
     m_spec.auto_stride(xstride, ystride, zstride, format, m_spec.nchannels,
                        m_spec.width, m_spec.height);
     const void* origdata = data;
+    if (format == TypeUnknown)
+        format = m_spec.format;
+
+    // PNG specifically dictates unassociated (un-"premultiplied") alpha.
+    // If we need to unassociate alpha, do it in float.
+    std::unique_ptr<float[]> unassoc_scratch;
+    size_t npixels = size_t(m_spec.width) * size_t(yend - ybegin);
+    size_t nvals   = npixels * size_t(m_spec.nchannels);
+    if (m_convert_alpha) {
+        unassoc_scratch.reset(new float[nvals]);
+        float* floatvals = unassoc_scratch.get();
+        // Contiguize and convert to float
+        OIIO::convert_image(m_spec.nchannels, m_spec.width, m_spec.height, 1,
+                            data, format, xstride, ystride, AutoStride,
+                            floatvals, TypeFloat, AutoStride, AutoStride,
+                            AutoStride);
+        // Deassociate alpha
+        deassociateAlpha(floatvals, npixels, m_spec.nchannels,
+                         m_spec.alpha_channel, m_gamma);
+        data    = floatvals;
+        format  = TypeFloat;
+        xstride = size_t(m_spec.nchannels) * sizeof(float);
+        ystride = xstride * size_t(m_spec.width);
+        zstride = ystride * size_t(m_spec.height);
+    }
+
     data = to_native_rectangle(m_spec.x, m_spec.x + m_spec.width, ybegin, yend,
                                z, z + 1, format, data, xstride, ystride,
-                               zstride, m_scratch);
-    size_t npixels = m_spec.width * (yend - ybegin);
-    size_t nvals   = npixels * m_spec.nchannels;
+                               zstride, m_scratch, m_dither, 0, ybegin, z);
     if (data == origdata && (m_convert_alpha || m_need_swap)) {
         m_scratch.assign((unsigned char*)data,
                          (unsigned char*)data + nvals * m_spec.format.size());
         data = m_scratch.data();
     }
 
-    // PNG specifically dictates unassociated (un-"premultiplied") alpha
-    if (m_convert_alpha) {
-        if (m_spec.format == TypeDesc::UINT16)
-            deassociateAlpha((unsigned short*)data, npixels, m_spec.nchannels,
-                             m_spec.alpha_channel, m_gamma);
-        else
-            deassociateAlpha((unsigned char*)data, npixels, m_spec.nchannels,
-                             m_spec.alpha_channel, m_gamma);
-    }
-
     // PNG is always big endian
     if (m_need_swap)
         swap_endian((unsigned short*)data, nvals);

testsuite/png/ref/out-libpng15.txt

@@ -27,3 +27,7 @@ exif.png             :   64 x   64, 3 channel, uint8 png
     SHA-1: 7CB41FEA50720B48BE0C145E1473982B23E9AB77
     channel list: R, G, B
     oiio:ColorSpace: "sRGB"
+smallalpha.png       :    1 x    1, 4 channel, uint8 png
+    Pixel (0, 0): 240 108 119 1 (0.94117653 0.42352945 0.4666667 0.003921569)
+Comparing "test16.png" and "ref/test16.png"
+PASS

testsuite/png/ref/out.txt

@@ -31,3 +31,7 @@ exif.png             :   64 x   64, 3 channel, uint8 png
     Exif:FocalLength: 45.7 (45.7 mm)
     Exif:WhiteBalance: 0 (auto)
     oiio:ColorSpace: "sRGB"
+smallalpha.png       :    1 x    1, 4 channel, uint8 png
+    Pixel (0, 0): 240 108 119 1 (0.94117653 0.42352945 0.4666667 0.003921569)
+Comparing "test16.png" and "ref/test16.png"
+PASS

testsuite/png/ref/test16.png

@@ -17,3 +17,14 @@
                      "--attrib Exif:FocalLength 45.7 " +
                      "-o exif.png")
 command += info_command ("exif.png", safematch=True)
+
+# regression test for 16 bit output bug
+command += oiiotool ("--pattern fill:topleft=1,0,0,1:topright=0,1,0,1:bottomleft=0,0,1,1:bottomright=1,1,1,1 16x16 4 -d uint16 -o test16.png")
+
+# Test high quality alpha deassociation using alpha value close to zero.
+# This example is inspired by Yafes on the Slack.
+command += oiiotool ("--pattern fill:color=0.00235,0.00106,0.00117,0.0025 1x1 4 -d uint8 -o smallalpha.png")
+command += oiiotool ("--no-autopremult --dumpdata smallalpha.png")
+
+outputs = [ "test16.png", "out.txt" ]
+