src/lib/image.cc

   1 /*
   2     Copyright (C) 2012-2021 Carl Hetherington <cth@carlh.net>
   3
   4     This file is part of DCP-o-matic.
   5
   6     DCP-o-matic is free software; you can redistribute it and/or modify
   7     it under the terms of the GNU General Public License as published by
   8     the Free Software Foundation; either version 2 of the License, or
   9     (at your option) any later version.
  10
  11     DCP-o-matic is distributed in the hope that it will be useful,
  12     but WITHOUT ANY WARRANTY; without even the implied warranty of
  13     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14     GNU General Public License for more details.
  15
  16     You should have received a copy of the GNU General Public License
  17     along with DCP-o-matic.  If not, see <http://www.gnu.org/licenses/>.
  18
  19 */
  20
  21
  22 /** @file src/image.cc
  23  *  @brief A class to describe a video image.
  24  */
  25
  26
  27 #include "compose.hpp"
  28 #include "dcpomatic_socket.h"
  29 #include "exceptions.h"
  30 #include "image.h"
  31 #include "rect.h"
  32 #include "timer.h"
  33 #include "util.h"
  34 #include <dcp/rgb_xyz.h>
  35 #include <dcp/transfer_function.h>
  36 extern "C" {
  37 #include <libavutil/frame.h>
  38 #include <libavutil/pixdesc.h>
  39 #include <libavutil/pixfmt.h>
  40 #include <libswscale/swscale.h>
  41 }
  42 #include <png.h>
  43 #if HAVE_VALGRIND_MEMCHECK_H
  44 #include <valgrind/memcheck.h>
  45 #endif
  46 #include <iostream>
  47
  48
  49 #include "i18n.h"
  50
  51
  52 using std::cerr;
  53 using std::cout;
  54 using std::list;
  55 using std::make_shared;
  56 using std::max;
  57 using std::min;
  58 using std::runtime_error;
  59 using std::shared_ptr;
  60 using std::string;
  61 using dcp::Size;
  62
  63
  64 /** The memory alignment, in bytes, used for each row of an image if aligment is requested */
  65 #define ALIGNMENT 64
  66
  67 /* U/V black value for 8-bit colour */
  68 static uint8_t const eight_bit_uv =     (1 << 7) - 1;
  69 /* U/V black value for 9-bit colour */
  70 static uint16_t const nine_bit_uv =     (1 << 8) - 1;
  71 /* U/V black value for 10-bit colour */
  72 static uint16_t const ten_bit_uv =      (1 << 9) - 1;
  73 /* U/V black value for 16-bit colour */
  74 static uint16_t const sixteen_bit_uv =  (1 << 15) - 1;
  75
  76
  77 int
  78 Image::vertical_factor (int n) const
  79 {
  80         if (n == 0) {
  81                 return 1;
  82         }
  83
  84         auto d = av_pix_fmt_desc_get(_pixel_format);
  85         if (!d) {
  86                 throw PixelFormatError ("line_factor()", _pixel_format);
  87         }
  88
  89         return lrintf(powf(2.0f, d->log2_chroma_h));
  90 }
  91
  92 int
  93 Image::horizontal_factor (int n) const
  94 {
  95         if (n == 0) {
  96                 return 1;
  97         }
  98
  99         auto d = av_pix_fmt_desc_get(_pixel_format);
 100         if (!d) {
 101                 throw PixelFormatError ("sample_size()", _pixel_format);
 102         }
 103
 104         return lrintf(powf(2.0f, d->log2_chroma_w));
 105 }
 106
 107 /** @param n Component index.
 108  *  @return Number of samples (i.e. pixels, unless sub-sampled) in each direction for this component.
 109  */
 110 dcp::Size
 111 Image::sample_size (int n) const
 112 {
 113         return dcp::Size (
 114                 lrint (ceil(static_cast<double>(size().width) / horizontal_factor(n))),
 115                 lrint (ceil(static_cast<double>(size().height) / vertical_factor(n)))
 116                 );
 117 }
 118
 119 /** @return Number of planes */
 120 int
 121 Image::planes () const
 122 {
 123         auto d = av_pix_fmt_desc_get(_pixel_format);
 124         if (!d) {
 125                 throw PixelFormatError ("planes()", _pixel_format);
 126         }
 127
 128         if (_pixel_format == AV_PIX_FMT_PAL8) {
 129                 return 2;
 130         }
 131
 132         if ((d->flags & AV_PIX_FMT_FLAG_PLANAR) == 0) {
 133                 return 1;
 134         }
 135
 136         return d->nb_components;
 137 }
 138
 139
 140 static
 141 int
 142 round_width_for_subsampling (int p, AVPixFmtDescriptor const * desc)
 143 {
 144         return p & ~ ((1 << desc->log2_chroma_w) - 1);
 145 }
 146
 147
 148 static
 149 int
 150 round_height_for_subsampling (int p, AVPixFmtDescriptor const * desc)
 151 {
 152         return p & ~ ((1 << desc->log2_chroma_h) - 1);
 153 }
 154
 155
 156 /** Crop this image, scale it to `inter_size' and then place it in a black frame of `out_size'.
 157  *  @param crop Amount to crop by.
 158  *  @param inter_size Size to scale the cropped image to.
 159  *  @param out_size Size of output frame; if this is larger than inter_size there will be black padding.
 160  *  @param yuv_to_rgb YUV to RGB transformation to use, if required.
 161  *  @param video_range Video range of the image.
 162  *  @param out_format Output pixel format.
 163  *  @param out_aligned true to make the output image aligned.
 164  *  @param out_video_range Video range to use for the output image.
 165  *  @param fast Try to be fast at the possible expense of quality; at present this means using
 166  *  fast bilinear rather than bicubic scaling.
 167  */
 168 shared_ptr<Image>
 169 Image::crop_scale_window (
 170         Crop crop,
 171         dcp::Size inter_size,
 172         dcp::Size out_size,
 173         dcp::YUVToRGB yuv_to_rgb,
 174         VideoRange video_range,
 175         AVPixelFormat out_format,
 176         VideoRange out_video_range,
 177         bool out_aligned,
 178         bool fast
 179         ) const
 180 {
 181         /* Empirical testing suggests that sws_scale() will crash if
 182            the input image is not aligned.
 183         */
 184         DCPOMATIC_ASSERT (aligned ());
 185
 186         DCPOMATIC_ASSERT (out_size.width >= inter_size.width);
 187         DCPOMATIC_ASSERT (out_size.height >= inter_size.height);
 188
 189         auto out = make_shared<Image>(out_format, out_size, out_aligned);
 190         out->make_black ();
 191
 192         auto in_desc = av_pix_fmt_desc_get (_pixel_format);
 193         if (!in_desc) {
 194                 throw PixelFormatError ("crop_scale_window()", _pixel_format);
 195         }
 196
 197         /* Round down so that we crop only the number of pixels that is straightforward
 198          * considering any subsampling.
 199          */
 200         Crop corrected_crop(
 201                 round_width_for_subsampling(crop.left, in_desc),
 202                 round_width_for_subsampling(crop.right, in_desc),
 203                 round_height_for_subsampling(crop.top, in_desc),
 204                 round_height_for_subsampling(crop.bottom, in_desc)
 205                 );
 206
 207         /* Also check that we aren't cropping more image than there actually is */
 208         if ((corrected_crop.left + corrected_crop.right) >= (size().width - 4)) {
 209                 corrected_crop.left = 0;
 210                 corrected_crop.right = size().width - 4;
 211         }
 212
 213         if ((corrected_crop.top + corrected_crop.bottom) >= (size().height - 4)) {
 214                 corrected_crop.top = 0;
 215                 corrected_crop.bottom = size().height - 4;
 216         }
 217
 218         /* Size of the image after any crop */
 219         auto const cropped_size = corrected_crop.apply (size());
 220
 221         /* Scale context for a scale from cropped_size to inter_size */
 222         auto scale_context = sws_getContext (
 223                         cropped_size.width, cropped_size.height, pixel_format(),
 224                         inter_size.width, inter_size.height, out_format,
 225                         fast ? SWS_FAST_BILINEAR : SWS_BICUBIC, 0, 0, 0
 226                 );
 227
 228         if (!scale_context) {
 229                 throw runtime_error (N_("Could not allocate SwsContext"));
 230         }
 231
 232         DCPOMATIC_ASSERT (yuv_to_rgb < dcp::YUVToRGB::COUNT);
 233         int const lut[static_cast<int>(dcp::YUVToRGB::COUNT)] = {
 234                 SWS_CS_ITU601,
 235                 SWS_CS_ITU709
 236         };
 237
 238         /* The 3rd parameter here is:
 239            0 -> source range MPEG (i.e. "video", 16-235)
 240            1 -> source range JPEG (i.e. "full", 0-255)
 241            And the 5th:
 242            0 -> destination range MPEG (i.e. "video", 16-235)
 243            1 -> destination range JPEG (i.e. "full", 0-255)
 244
 245            But remember: sws_setColorspaceDetails ignores these
 246            parameters unless the both source and destination images
 247            are isYUV or isGray.  (If either is not, it uses video range).
 248         */
 249         sws_setColorspaceDetails (
 250                 scale_context,
 251                 sws_getCoefficients (lut[static_cast<int>(yuv_to_rgb)]), video_range == VideoRange::VIDEO ? 0 : 1,
 252                 sws_getCoefficients (lut[static_cast<int>(yuv_to_rgb)]), out_video_range == VideoRange::VIDEO ? 0 : 1,
 253                 0, 1 << 16, 1 << 16
 254                 );
 255
 256         /* Prepare input data pointers with crop */
 257         uint8_t* scale_in_data[planes()];
 258         for (int c = 0; c < planes(); ++c) {
 259                 int const x = lrintf(bytes_per_pixel(c) * corrected_crop.left);
 260                 scale_in_data[c] = data()[c] + x + stride()[c] * (corrected_crop.top / vertical_factor(c));
 261         }
 262
 263         auto out_desc = av_pix_fmt_desc_get (out_format);
 264         if (!out_desc) {
 265                 throw PixelFormatError ("crop_scale_window()", out_format);
 266         }
 267
 268         /* Corner of the image within out_size */
 269         Position<int> const corner (
 270                 round_width_for_subsampling((out_size.width - inter_size.width) / 2, out_desc),
 271                 round_height_for_subsampling((out_size.height - inter_size.height) / 2, out_desc)
 272                 );
 273
 274         uint8_t* scale_out_data[out->planes()];
 275         for (int c = 0; c < out->planes(); ++c) {
 276                 int const x = lrintf(out->bytes_per_pixel(c) * corner.x);
 277                 scale_out_data[c] = out->data()[c] + x + out->stride()[c] * (corner.y / out->vertical_factor(c));
 278         }
 279
 280         sws_scale (
 281                 scale_context,
 282                 scale_in_data, stride(),
 283                 0, cropped_size.height,
 284                 scale_out_data, out->stride()
 285                 );
 286
 287         sws_freeContext (scale_context);
 288
 289         if (corrected_crop != Crop() && cropped_size == inter_size) {
 290                 /* We are cropping without any scaling or pixel format conversion, so FFmpeg may have left some
 291                    data behind in our image.  Clear it out.  It may get to the point where we should just stop
 292                    trying to be clever with cropping.
 293                 */
 294                 out->make_part_black (corner.x + cropped_size.width, out_size.width - cropped_size.width);
 295         }
 296
 297         return out;
 298 }
 299
 300 shared_ptr<Image>
 301 Image::convert_pixel_format (dcp::YUVToRGB yuv_to_rgb, AVPixelFormat out_format, bool out_aligned, bool fast) const
 302 {
 303         return scale(size(), yuv_to_rgb, out_format, out_aligned, fast);
 304 }
 305
 306 /** @param out_size Size to scale to.
 307  *  @param yuv_to_rgb YUVToRGB transform transform to use, if required.
 308  *  @param out_format Output pixel format.
 309  *  @param out_aligned true to make an aligned output image.
 310  *  @param fast Try to be fast at the possible expense of quality; at present this means using
 311  *  fast bilinear rather than bicubic scaling.
 312  */
 313 shared_ptr<Image>
 314 Image::scale (dcp::Size out_size, dcp::YUVToRGB yuv_to_rgb, AVPixelFormat out_format, bool out_aligned, bool fast) const
 315 {
 316         /* Empirical testing suggests that sws_scale() will crash if
 317            the input image is not aligned.
 318         */
 319         DCPOMATIC_ASSERT (aligned ());
 320
 321         auto scaled = make_shared<Image>(out_format, out_size, out_aligned);
 322         auto scale_context = sws_getContext (
 323                 size().width, size().height, pixel_format(),
 324                 out_size.width, out_size.height, out_format,
 325                 (fast ? SWS_FAST_BILINEAR : SWS_BICUBIC) | SWS_ACCURATE_RND, 0, 0, 0
 326                 );
 327
 328         DCPOMATIC_ASSERT (yuv_to_rgb < dcp::YUVToRGB::COUNT);
 329         int const lut[static_cast<int>(dcp::YUVToRGB::COUNT)] = {
 330                 SWS_CS_ITU601,
 331                 SWS_CS_ITU709
 332         };
 333
 334         /* The 3rd parameter here is:
 335            0 -> source range MPEG (i.e. "video", 16-235)
 336            1 -> source range JPEG (i.e. "full", 0-255)
 337            And the 5th:
 338            0 -> destination range MPEG (i.e. "video", 16-235)
 339            1 -> destination range JPEG (i.e. "full", 0-255)
 340
 341            But remember: sws_setColorspaceDetails ignores these
 342            parameters unless the corresponding image isYUV or isGray.
 343            (If it's neither, it uses video range).
 344         */
 345         sws_setColorspaceDetails (
 346                 scale_context,
 347                 sws_getCoefficients (lut[static_cast<int>(yuv_to_rgb)]), 0,
 348                 sws_getCoefficients (lut[static_cast<int>(yuv_to_rgb)]), 0,
 349                 0, 1 << 16, 1 << 16
 350                 );
 351
 352         sws_scale (
 353                 scale_context,
 354                 data(), stride(),
 355                 0, size().height,
 356                 scaled->data(), scaled->stride()
 357                 );
 358
 359         sws_freeContext (scale_context);
 360
 361         return scaled;
 362 }
 363
 364 /** Blacken a YUV image whose bits per pixel is rounded up to 16 */
 365 void
 366 Image::yuv_16_black (uint16_t v, bool alpha)
 367 {
 368         memset (data()[0], 0, sample_size(0).height * stride()[0]);
 369         for (int i = 1; i < 3; ++i) {
 370                 auto p = reinterpret_cast<int16_t*> (data()[i]);
 371                 int const lines = sample_size(i).height;
 372                 for (int y = 0; y < lines; ++y) {
 373                         /* We divide by 2 here because we are writing 2 bytes at a time */
 374                         for (int x = 0; x < line_size()[i] / 2; ++x) {
 375                                 p[x] = v;
 376                         }
 377                         p += stride()[i] / 2;
 378                 }
 379         }
 380
 381         if (alpha) {
 382                 memset (data()[3], 0, sample_size(3).height * stride()[3]);
 383         }
 384 }
 385
 386 uint16_t
 387 Image::swap_16 (uint16_t v)
 388 {
 389         return ((v >> 8) & 0xff) | ((v & 0xff) << 8);
 390 }
 391
 392 void
 393 Image::make_part_black (int const start, int const width)
 394 {
 395         auto y_part = [&]() {
 396                 int const bpp = bytes_per_pixel(0);
 397                 int const h = sample_size(0).height;
 398                 int const s = stride()[0];
 399                 auto p = data()[0];
 400                 for (int y = 0; y < h; ++y) {
 401                         memset (p + start * bpp, 0, width * bpp);
 402                         p += s;
 403                 }
 404         };
 405
 406         switch (_pixel_format) {
 407         case AV_PIX_FMT_RGB24:
 408         case AV_PIX_FMT_ARGB:
 409         case AV_PIX_FMT_RGBA:
 410         case AV_PIX_FMT_ABGR:
 411         case AV_PIX_FMT_BGRA:
 412         case AV_PIX_FMT_RGB555LE:
 413         case AV_PIX_FMT_RGB48LE:
 414         case AV_PIX_FMT_RGB48BE:
 415         case AV_PIX_FMT_XYZ12LE:
 416         {
 417                 int const h = sample_size(0).height;
 418                 int const bpp = bytes_per_pixel(0);
 419                 int const s = stride()[0];
 420                 uint8_t* p = data()[0];
 421                 for (int y = 0; y < h; y++) {
 422                         memset (p + start * bpp, 0, width * bpp);
 423                         p += s;
 424                 }
 425                 break;
 426         }
 427         case AV_PIX_FMT_YUV420P:
 428         {
 429                 y_part ();
 430                 for (int i = 1; i < 3; ++i) {
 431                         auto p = data()[i];
 432                         int const h = sample_size(i).height;
 433                         for (int y = 0; y < h; ++y) {
 434                                 for (int x = start / 2; x < (start + width) / 2; ++x) {
 435                                         p[x] = eight_bit_uv;
 436                                 }
 437                                 p += stride()[i];
 438                         }
 439                 }
 440                 break;
 441         }
 442         case AV_PIX_FMT_YUV422P10LE:
 443         {
 444                 y_part ();
 445                 for (int i = 1; i < 3; ++i) {
 446                         auto p = reinterpret_cast<int16_t*>(data()[i]);
 447                         int const h = sample_size(i).height;
 448                         for (int y = 0; y < h; ++y) {
 449                                 for (int x = start / 2; x < (start + width) / 2; ++x) {
 450                                         p[x] = ten_bit_uv;
 451                                 }
 452                                 p += stride()[i] / 2;
 453                         }
 454                 }
 455                 break;
 456         }
 457         default:
 458                 throw PixelFormatError ("make_part_black()", _pixel_format);
 459         }
 460 }
 461
 462 void
 463 Image::make_black ()
 464 {
 465         switch (_pixel_format) {
 466         case AV_PIX_FMT_YUV420P:
 467         case AV_PIX_FMT_YUV422P:
 468         case AV_PIX_FMT_YUV444P:
 469         case AV_PIX_FMT_YUV411P:
 470                 memset (data()[0], 0, sample_size(0).height * stride()[0]);
 471                 memset (data()[1], eight_bit_uv, sample_size(1).height * stride()[1]);
 472                 memset (data()[2], eight_bit_uv, sample_size(2).height * stride()[2]);
 473                 break;
 474
 475         case AV_PIX_FMT_YUVJ420P:
 476         case AV_PIX_FMT_YUVJ422P:
 477         case AV_PIX_FMT_YUVJ444P:
 478                 memset (data()[0], 0, sample_size(0).height * stride()[0]);
 479                 memset (data()[1], eight_bit_uv + 1, sample_size(1).height * stride()[1]);
 480                 memset (data()[2], eight_bit_uv + 1, sample_size(2).height * stride()[2]);
 481                 break;
 482
 483         case AV_PIX_FMT_YUV422P9LE:
 484         case AV_PIX_FMT_YUV444P9LE:
 485                 yuv_16_black (nine_bit_uv, false);
 486                 break;
 487
 488         case AV_PIX_FMT_YUV422P9BE:
 489         case AV_PIX_FMT_YUV444P9BE:
 490                 yuv_16_black (swap_16 (nine_bit_uv), false);
 491                 break;
 492
 493         case AV_PIX_FMT_YUV422P10LE:
 494         case AV_PIX_FMT_YUV444P10LE:
 495                 yuv_16_black (ten_bit_uv, false);
 496                 break;
 497
 498         case AV_PIX_FMT_YUV422P16LE:
 499         case AV_PIX_FMT_YUV444P16LE:
 500                 yuv_16_black (sixteen_bit_uv, false);
 501                 break;
 502
 503         case AV_PIX_FMT_YUV444P10BE:
 504         case AV_PIX_FMT_YUV422P10BE:
 505                 yuv_16_black (swap_16 (ten_bit_uv), false);
 506                 break;
 507
 508         case AV_PIX_FMT_YUVA420P9BE:
 509         case AV_PIX_FMT_YUVA422P9BE:
 510         case AV_PIX_FMT_YUVA444P9BE:
 511                 yuv_16_black (swap_16 (nine_bit_uv), true);
 512                 break;
 513
 514         case AV_PIX_FMT_YUVA420P9LE:
 515         case AV_PIX_FMT_YUVA422P9LE:
 516         case AV_PIX_FMT_YUVA444P9LE:
 517                 yuv_16_black (nine_bit_uv, true);
 518                 break;
 519
 520         case AV_PIX_FMT_YUVA420P10BE:
 521         case AV_PIX_FMT_YUVA422P10BE:
 522         case AV_PIX_FMT_YUVA444P10BE:
 523                 yuv_16_black (swap_16 (ten_bit_uv), true);
 524                 break;
 525
 526         case AV_PIX_FMT_YUVA420P10LE:
 527         case AV_PIX_FMT_YUVA422P10LE:
 528         case AV_PIX_FMT_YUVA444P10LE:
 529                 yuv_16_black (ten_bit_uv, true);
 530                 break;
 531
 532         case AV_PIX_FMT_YUVA420P16BE:
 533         case AV_PIX_FMT_YUVA422P16BE:
 534         case AV_PIX_FMT_YUVA444P16BE:
 535                 yuv_16_black (swap_16 (sixteen_bit_uv), true);
 536                 break;
 537
 538         case AV_PIX_FMT_YUVA420P16LE:
 539         case AV_PIX_FMT_YUVA422P16LE:
 540         case AV_PIX_FMT_YUVA444P16LE:
 541                 yuv_16_black (sixteen_bit_uv, true);
 542                 break;
 543
 544         case AV_PIX_FMT_RGB24:
 545         case AV_PIX_FMT_ARGB:
 546         case AV_PIX_FMT_RGBA:
 547         case AV_PIX_FMT_ABGR:
 548         case AV_PIX_FMT_BGRA:
 549         case AV_PIX_FMT_RGB555LE:
 550         case AV_PIX_FMT_RGB48LE:
 551         case AV_PIX_FMT_RGB48BE:
 552         case AV_PIX_FMT_XYZ12LE:
 553                 memset (data()[0], 0, sample_size(0).height * stride()[0]);
 554                 break;
 555
 556         case AV_PIX_FMT_UYVY422:
 557         {
 558                 int const Y = sample_size(0).height;
 559                 int const X = line_size()[0];
 560                 uint8_t* p = data()[0];
 561                 for (int y = 0; y < Y; ++y) {
 562                         for (int x = 0; x < X / 4; ++x) {
 563                                 *p++ = eight_bit_uv; // Cb
 564                                 *p++ = 0;            // Y0
 565                                 *p++ = eight_bit_uv; // Cr
 566                                 *p++ = 0;            // Y1
 567                         }
 568                 }
 569                 break;
 570         }
 571
 572         default:
 573                 throw PixelFormatError ("make_black()", _pixel_format);
 574         }
 575 }
 576
 577 void
 578 Image::make_transparent ()
 579 {
 580         if (_pixel_format != AV_PIX_FMT_BGRA && _pixel_format != AV_PIX_FMT_RGBA) {
 581                 throw PixelFormatError ("make_transparent()", _pixel_format);
 582         }
 583
 584         memset (data()[0], 0, sample_size(0).height * stride()[0]);
 585 }
 586
 587 void
 588 Image::alpha_blend (shared_ptr<const Image> other, Position<int> position)
 589 {
 590         /* We're blending RGBA or BGRA images */
 591         DCPOMATIC_ASSERT (other->pixel_format() == AV_PIX_FMT_BGRA || other->pixel_format() == AV_PIX_FMT_RGBA);
 592         int const blue = other->pixel_format() == AV_PIX_FMT_BGRA ? 0 : 2;
 593         int const red = other->pixel_format() == AV_PIX_FMT_BGRA ? 2 : 0;
 594
 595         int const other_bpp = 4;
 596
 597         int start_tx = position.x;
 598         int start_ox = 0;
 599
 600         if (start_tx < 0) {
 601                 start_ox = -start_tx;
 602                 start_tx = 0;
 603         }
 604
 605         int start_ty = position.y;
 606         int start_oy = 0;
 607
 608         if (start_ty < 0) {
 609                 start_oy = -start_ty;
 610                 start_ty = 0;
 611         }
 612
 613         switch (_pixel_format) {
 614         case AV_PIX_FMT_RGB24:
 615         {
 616                 /* Going onto RGB24.  First byte is red, second green, third blue */
 617                 int const this_bpp = 3;
 618                 for (int ty = start_ty, oy = start_oy; ty < size().height && oy < other->size().height; ++ty, ++oy) {
 619                         uint8_t* tp = data()[0] + ty * stride()[0] + start_tx * this_bpp;
 620                         uint8_t* op = other->data()[0] + oy * other->stride()[0];
 621                         for (int tx = start_tx, ox = start_ox; tx < size().width && ox < other->size().width; ++tx, ++ox) {
 622                                 float const alpha = float (op[3]) / 255;
 623                                 tp[0] = op[red] * alpha + tp[0] * (1 - alpha);
 624                                 tp[1] = op[1] * alpha + tp[1] * (1 - alpha);
 625                                 tp[2] = op[blue] * alpha + tp[2] * (1 - alpha);
 626
 627                                 tp += this_bpp;
 628                                 op += other_bpp;
 629                         }
 630                 }
 631                 break;
 632         }
 633         case AV_PIX_FMT_BGRA:
 634         {
 635                 int const this_bpp = 4;
 636                 for (int ty = start_ty, oy = start_oy; ty < size().height && oy < other->size().height; ++ty, ++oy) {
 637                         uint8_t* tp = data()[0] + ty * stride()[0] + start_tx * this_bpp;
 638                         uint8_t* op = other->data()[0] + oy * other->stride()[0];
 639                         for (int tx = start_tx, ox = start_ox; tx < size().width && ox < other->size().width; ++tx, ++ox) {
 640                                 float const alpha = float (op[3]) / 255;
 641                                 tp[0] = op[blue] * alpha + tp[0] * (1 - alpha);
 642                                 tp[1] = op[1] * alpha + tp[1] * (1 - alpha);
 643                                 tp[2] = op[red] * alpha + tp[2] * (1 - alpha);
 644                                 tp[3] = op[3] * alpha + tp[3] * (1 - alpha);
 645
 646                                 tp += this_bpp;
 647                                 op += other_bpp;
 648                         }
 649                 }
 650                 break;
 651         }
 652         case AV_PIX_FMT_RGBA:
 653         {
 654                 int const this_bpp = 4;
 655                 for (int ty = start_ty, oy = start_oy; ty < size().height && oy < other->size().height; ++ty, ++oy) {
 656                         uint8_t* tp = data()[0] + ty * stride()[0] + start_tx * this_bpp;
 657                         uint8_t* op = other->data()[0] + oy * other->stride()[0];
 658                         for (int tx = start_tx, ox = start_ox; tx < size().width && ox < other->size().width; ++tx, ++ox) {
 659                                 float const alpha = float (op[3]) / 255;
 660                                 tp[0] = op[red] * alpha + tp[0] * (1 - alpha);
 661                                 tp[1] = op[1] * alpha + tp[1] * (1 - alpha);
 662                                 tp[2] = op[blue] * alpha + tp[2] * (1 - alpha);
 663                                 tp[3] = op[3] * alpha + tp[3] * (1 - alpha);
 664
 665                                 tp += this_bpp;
 666                                 op += other_bpp;
 667                         }
 668                 }
 669                 break;
 670         }
 671         case AV_PIX_FMT_RGB48LE:
 672         {
 673                 int const this_bpp = 6;
 674                 for (int ty = start_ty, oy = start_oy; ty < size().height && oy < other->size().height; ++ty, ++oy) {
 675                         uint8_t* tp = data()[0] + ty * stride()[0] + start_tx * this_bpp;
 676                         uint8_t* op = other->data()[0] + oy * other->stride()[0];
 677                         for (int tx = start_tx, ox = start_ox; tx < size().width && ox < other->size().width; ++tx, ++ox) {
 678                                 float const alpha = float (op[3]) / 255;
 679                                 /* Blend high bytes */
 680                                 tp[1] = op[red] * alpha + tp[1] * (1 - alpha);
 681                                 tp[3] = op[1] * alpha + tp[3] * (1 - alpha);
 682                                 tp[5] = op[blue] * alpha + tp[5] * (1 - alpha);
 683
 684                                 tp += this_bpp;
 685                                 op += other_bpp;
 686                         }
 687                 }
 688                 break;
 689         }
 690         case AV_PIX_FMT_XYZ12LE:
 691         {
 692                 auto conv = dcp::ColourConversion::srgb_to_xyz();
 693                 double fast_matrix[9];
 694                 dcp::combined_rgb_to_xyz (conv, fast_matrix);
 695                 double const * lut_in = conv.in()->lut (8, false);
 696                 double const * lut_out = conv.out()->lut (16, true);
 697                 int const this_bpp = 6;
 698                 for (int ty = start_ty, oy = start_oy; ty < size().height && oy < other->size().height; ++ty, ++oy) {
 699                         uint16_t* tp = reinterpret_cast<uint16_t*> (data()[0] + ty * stride()[0] + start_tx * this_bpp);
 700                         uint8_t* op = other->data()[0] + oy * other->stride()[0];
 701                         for (int tx = start_tx, ox = start_ox; tx < size().width && ox < other->size().width; ++tx, ++ox) {
 702                                 float const alpha = float (op[3]) / 255;
 703
 704                                 /* Convert sRGB to XYZ; op is BGRA.  First, input gamma LUT */
 705                                 double const r = lut_in[op[red]];
 706                                 double const g = lut_in[op[1]];
 707                                 double const b = lut_in[op[blue]];
 708
 709                                 /* RGB to XYZ, including Bradford transform and DCI companding */
 710                                 double const x = max (0.0, min (65535.0, r * fast_matrix[0] + g * fast_matrix[1] + b * fast_matrix[2]));
 711                                 double const y = max (0.0, min (65535.0, r * fast_matrix[3] + g * fast_matrix[4] + b * fast_matrix[5]));
 712                                 double const z = max (0.0, min (65535.0, r * fast_matrix[6] + g * fast_matrix[7] + b * fast_matrix[8]));
 713
 714                                 /* Out gamma LUT and blend */
 715                                 tp[0] = lrint(lut_out[lrint(x)] * 65535) * alpha + tp[0] * (1 - alpha);
 716                                 tp[1] = lrint(lut_out[lrint(y)] * 65535) * alpha + tp[1] * (1 - alpha);
 717                                 tp[2] = lrint(lut_out[lrint(z)] * 65535) * alpha + tp[2] * (1 - alpha);
 718
 719                                 tp += this_bpp / 2;
 720                                 op += other_bpp;
 721                         }
 722                 }
 723                 break;
 724         }
 725         case AV_PIX_FMT_YUV420P:
 726         {
 727                 auto yuv = other->convert_pixel_format (dcp::YUVToRGB::REC709, _pixel_format, false, false);
 728                 dcp::Size const ts = size();
 729                 dcp::Size const os = yuv->size();
 730                 for (int ty = start_ty, oy = start_oy; ty < ts.height && oy < os.height; ++ty, ++oy) {
 731                         int const hty = ty / 2;
 732                         int const hoy = oy / 2;
 733                         uint8_t* tY = data()[0] + (ty * stride()[0]) + start_tx;
 734                         uint8_t* tU = data()[1] + (hty * stride()[1]) + start_tx / 2;
 735                         uint8_t* tV = data()[2] + (hty * stride()[2]) + start_tx / 2;
 736                         uint8_t* oY = yuv->data()[0] + (oy * yuv->stride()[0]) + start_ox;
 737                         uint8_t* oU = yuv->data()[1] + (hoy * yuv->stride()[1]) + start_ox / 2;
 738                         uint8_t* oV = yuv->data()[2] + (hoy * yuv->stride()[2]) + start_ox / 2;
 739                         uint8_t* alpha = other->data()[0] + (oy * other->stride()[0]) + start_ox * 4;
 740                         for (int tx = start_tx, ox = start_ox; tx < ts.width && ox < os.width; ++tx, ++ox) {
 741                                 float const a = float(alpha[3]) / 255;
 742                                 *tY = *oY * a + *tY * (1 - a);
 743                                 *tU = *oU * a + *tU * (1 - a);
 744                                 *tV = *oV * a + *tV * (1 - a);
 745                                 ++tY;
 746                                 ++oY;
 747                                 if (tx % 2) {
 748                                         ++tU;
 749                                         ++tV;
 750                                 }
 751                                 if (ox % 2) {
 752                                         ++oU;
 753                                         ++oV;
 754                                 }
 755                                 alpha += 4;
 756                         }
 757                 }
 758                 break;
 759         }
 760         case AV_PIX_FMT_YUV420P10:
 761         {
 762                 auto yuv = other->convert_pixel_format (dcp::YUVToRGB::REC709, _pixel_format, false, false);
 763                 dcp::Size const ts = size();
 764                 dcp::Size const os = yuv->size();
 765                 for (int ty = start_ty, oy = start_oy; ty < ts.height && oy < os.height; ++ty, ++oy) {
 766                         int const hty = ty / 2;
 767                         int const hoy = oy / 2;
 768                         uint16_t* tY = ((uint16_t *) (data()[0] + (ty * stride()[0]))) + start_tx;
 769                         uint16_t* tU = ((uint16_t *) (data()[1] + (hty * stride()[1]))) + start_tx / 2;
 770                         uint16_t* tV = ((uint16_t *) (data()[2] + (hty * stride()[2]))) + start_tx / 2;
 771                         uint16_t* oY = ((uint16_t *) (yuv->data()[0] + (oy * yuv->stride()[0]))) + start_ox;
 772                         uint16_t* oU = ((uint16_t *) (yuv->data()[1] + (hoy * yuv->stride()[1]))) + start_ox / 2;
 773                         uint16_t* oV = ((uint16_t *) (yuv->data()[2] + (hoy * yuv->stride()[2]))) + start_ox / 2;
 774                         uint8_t* alpha = other->data()[0] + (oy * other->stride()[0]) + start_ox * 4;
 775                         for (int tx = start_tx, ox = start_ox; tx < ts.width && ox < os.width; ++tx, ++ox) {
 776                                 float const a = float(alpha[3]) / 255;
 777                                 *tY = *oY * a + *tY * (1 - a);
 778                                 *tU = *oU * a + *tU * (1 - a);
 779                                 *tV = *oV * a + *tV * (1 - a);
 780                                 ++tY;
 781                                 ++oY;
 782                                 if (tx % 2) {
 783                                         ++tU;
 784                                         ++tV;
 785                                 }
 786                                 if (ox % 2) {
 787                                         ++oU;
 788                                         ++oV;
 789                                 }
 790                                 alpha += 4;
 791                         }
 792                 }
 793                 break;
 794         }
 795         case AV_PIX_FMT_YUV422P10LE:
 796         {
 797                 auto yuv = other->convert_pixel_format (dcp::YUVToRGB::REC709, _pixel_format, false, false);
 798                 dcp::Size const ts = size();
 799                 dcp::Size const os = yuv->size();
 800                 for (int ty = start_ty, oy = start_oy; ty < ts.height && oy < os.height; ++ty, ++oy) {
 801                         uint16_t* tY = ((uint16_t *) (data()[0] + (ty * stride()[0]))) + start_tx;
 802                         uint16_t* tU = ((uint16_t *) (data()[1] + (ty * stride()[1]))) + start_tx / 2;
 803                         uint16_t* tV = ((uint16_t *) (data()[2] + (ty * stride()[2]))) + start_tx / 2;
 804                         uint16_t* oY = ((uint16_t *) (yuv->data()[0] + (oy * yuv->stride()[0]))) + start_ox;
 805                         uint16_t* oU = ((uint16_t *) (yuv->data()[1] + (oy * yuv->stride()[1]))) + start_ox / 2;
 806                         uint16_t* oV = ((uint16_t *) (yuv->data()[2] + (oy * yuv->stride()[2]))) + start_ox / 2;
 807                         uint8_t* alpha = other->data()[0] + (oy * other->stride()[0]) + start_ox * 4;
 808                         for (int tx = start_tx, ox = start_ox; tx < ts.width && ox < os.width; ++tx, ++ox) {
 809                                 float const a = float(alpha[3]) / 255;
 810                                 *tY = *oY * a + *tY * (1 - a);
 811                                 *tU = *oU * a + *tU * (1 - a);
 812                                 *tV = *oV * a + *tV * (1 - a);
 813                                 ++tY;
 814                                 ++oY;
 815                                 if (tx % 2) {
 816                                         ++tU;
 817                                         ++tV;
 818                                 }
 819                                 if (ox % 2) {
 820                                         ++oU;
 821                                         ++oV;
 822                                 }
 823                                 alpha += 4;
 824                         }
 825                 }
 826                 break;
 827         }
 828         default:
 829                 throw PixelFormatError ("alpha_blend()", _pixel_format);
 830         }
 831 }
 832
 833 void
 834 Image::copy (shared_ptr<const Image> other, Position<int> position)
 835 {
 836         /* Only implemented for RGB24 onto RGB24 so far */
 837         DCPOMATIC_ASSERT (_pixel_format == AV_PIX_FMT_RGB24 && other->pixel_format() == AV_PIX_FMT_RGB24);
 838         DCPOMATIC_ASSERT (position.x >= 0 && position.y >= 0);
 839
 840         int const N = min (position.x + other->size().width, size().width) - position.x;
 841         for (int ty = position.y, oy = 0; ty < size().height && oy < other->size().height; ++ty, ++oy) {
 842                 uint8_t * const tp = data()[0] + ty * stride()[0] + position.x * 3;
 843                 uint8_t * const op = other->data()[0] + oy * other->stride()[0];
 844                 memcpy (tp, op, N * 3);
 845         }
 846 }
 847
 848 void
 849 Image::read_from_socket (shared_ptr<Socket> socket)
 850 {
 851         for (int i = 0; i < planes(); ++i) {
 852                 uint8_t* p = data()[i];
 853                 int const lines = sample_size(i).height;
 854                 for (int y = 0; y < lines; ++y) {
 855                         socket->read (p, line_size()[i]);
 856                         p += stride()[i];
 857                 }
 858         }
 859 }
 860
 861 void
 862 Image::write_to_socket (shared_ptr<Socket> socket) const
 863 {
 864         for (int i = 0; i < planes(); ++i) {
 865                 uint8_t* p = data()[i];
 866                 int const lines = sample_size(i).height;
 867                 for (int y = 0; y < lines; ++y) {
 868                         socket->write (p, line_size()[i]);
 869                         p += stride()[i];
 870                 }
 871         }
 872 }
 873
 874 float
 875 Image::bytes_per_pixel (int c) const
 876 {
 877         auto d = av_pix_fmt_desc_get(_pixel_format);
 878         if (!d) {
 879                 throw PixelFormatError ("bytes_per_pixel()", _pixel_format);
 880         }
 881
 882         if (c >= planes()) {
 883                 return 0;
 884         }
 885
 886         float bpp[4] = { 0, 0, 0, 0 };
 887
 888 #ifdef DCPOMATIC_HAVE_AVCOMPONENTDESCRIPTOR_DEPTH_MINUS1
 889         bpp[0] = floor ((d->comp[0].depth_minus1 + 8) / 8);
 890         if (d->nb_components > 1) {
 891                 bpp[1] = floor ((d->comp[1].depth_minus1 + 8) / 8) / pow (2.0f, d->log2_chroma_w);
 892         }
 893         if (d->nb_components > 2) {
 894                 bpp[2] = floor ((d->comp[2].depth_minus1 + 8) / 8) / pow (2.0f, d->log2_chroma_w);
 895         }
 896         if (d->nb_components > 3) {
 897                 bpp[3] = floor ((d->comp[3].depth_minus1 + 8) / 8) / pow (2.0f, d->log2_chroma_w);
 898         }
 899 #else
 900         bpp[0] = floor ((d->comp[0].depth + 7) / 8);
 901         if (d->nb_components > 1) {
 902                 bpp[1] = floor ((d->comp[1].depth + 7) / 8) / pow (2.0f, d->log2_chroma_w);
 903         }
 904         if (d->nb_components > 2) {
 905                 bpp[2] = floor ((d->comp[2].depth + 7) / 8) / pow (2.0f, d->log2_chroma_w);
 906         }
 907         if (d->nb_components > 3) {
 908                 bpp[3] = floor ((d->comp[3].depth + 7) / 8) / pow (2.0f, d->log2_chroma_w);
 909         }
 910 #endif
 911
 912         if ((d->flags & AV_PIX_FMT_FLAG_PLANAR) == 0) {
 913                 /* Not planar; sum them up */
 914                 return bpp[0] + bpp[1] + bpp[2] + bpp[3];
 915         }
 916
 917         return bpp[c];
 918 }
 919
 920 /** Construct a Image of a given size and format, allocating memory
 921  *  as required.
 922  *
 923  *  @param p Pixel format.
 924  *  @param s Size in pixels.
 925  *  @param aligned true to make each row of this image aligned to a ALIGNMENT-byte boundary.
 926  */
 927 Image::Image (AVPixelFormat p, dcp::Size s, bool aligned)
 928         : _size (s)
 929         , _pixel_format (p)
 930         , _aligned (aligned)
 931 {
 932         allocate ();
 933 }
 934
 935 void
 936 Image::allocate ()
 937 {
 938         _data = (uint8_t **) wrapped_av_malloc (4 * sizeof (uint8_t *));
 939         _data[0] = _data[1] = _data[2] = _data[3] = 0;
 940
 941         _line_size = (int *) wrapped_av_malloc (4 * sizeof (int));
 942         _line_size[0] = _line_size[1] = _line_size[2] = _line_size[3] = 0;
 943
 944         _stride = (int *) wrapped_av_malloc (4 * sizeof (int));
 945         _stride[0] = _stride[1] = _stride[2] = _stride[3] = 0;
 946
 947         for (int i = 0; i < planes(); ++i) {
 948                 _line_size[i] = ceil (_size.width * bytes_per_pixel(i));
 949                 _stride[i] = stride_round_up (i, _line_size, _aligned ? ALIGNMENT : 1);
 950
 951                 /* The assembler function ff_rgb24ToY_avx (in libswscale/x86/input.asm)
 952                    uses a 16-byte fetch to read three bytes (R/G/B) of image data.
 953                    Hence on the last pixel of the last line it reads over the end of
 954                    the actual data by 1 byte.  If the width of an image is a multiple
 955                    of the stride alignment there will be no padding at the end of image lines.
 956                    OS X crashes on this illegal read, though other operating systems don't
 957                    seem to mind.  The nasty + 1 in this malloc makes sure there is always a byte
 958                    for that instruction to read safely.
 959
 960                    Further to the above, valgrind is now telling me that ff_rgb24ToY_ssse3
 961                    over-reads by more then _avx.  I can't follow the code to work out how much,
 962                    so I'll just over-allocate by ALIGNMENT bytes and have done with it.  Empirical
 963                    testing suggests that it works.
 964
 965                    In addition to these concerns, we may read/write as much as a whole extra line
 966                    at the end of each plane in cases where we are messing with offsets in order to
 967                    do pad or crop.  To solve this we over-allocate by an extra _stride[i] bytes.
 968
 969                    As an example: we may write to images starting at an offset so we get some padding.
 970                    Hence we want to write in the following pattern:
 971
 972                    block start   write start                                  line end
 973                    |..(padding)..|<------line-size------------->|..(padding)..|
 974                    |..(padding)..|<------line-size------------->|..(padding)..|
 975                    |..(padding)..|<------line-size------------->|..(padding)..|
 976
 977                    where line-size is of the smaller (inter_size) image and the full padded line length is that of
 978                    out_size.  To get things to work we have to tell FFmpeg that the stride is that of out_size.
 979                    However some parts of FFmpeg (notably rgb48Toxyz12 in swscale.c) process data for the full
 980                    specified *stride*.  This does not matter until we get to the last line:
 981
 982                    block start   write start                                  line end
 983                    |..(padding)..|<------line-size------------->|XXXwrittenXXX|
 984                    |XXXwrittenXXX|<------line-size------------->|XXXwrittenXXX|
 985                    |XXXwrittenXXX|<------line-size------------->|XXXwrittenXXXXXXwrittenXXX
 986                                                                                ^^^^ out of bounds
 987                 */
 988                 _data[i] = (uint8_t *) wrapped_av_malloc (_stride[i] * (sample_size(i).height + 1) + ALIGNMENT);
 989 #if HAVE_VALGRIND_MEMCHECK_H
 990                 /* The data between the end of the line size and the stride is undefined but processed by
 991                    libswscale, causing lots of valgrind errors.  Mark it all defined to quell these errors.
 992                 */
 993                 VALGRIND_MAKE_MEM_DEFINED (_data[i], _stride[i] * (sample_size(i).height + 1) + ALIGNMENT);
 994 #endif
 995         }
 996 }
 997
 998 Image::Image (Image const & other)
 999         : std::enable_shared_from_this<Image>(other)
1000         , _size (other._size)
1001         , _pixel_format (other._pixel_format)
1002         , _aligned (other._aligned)
1003 {
1004         allocate ();
1005
1006         for (int i = 0; i < planes(); ++i) {
1007                 uint8_t* p = _data[i];
1008                 uint8_t* q = other._data[i];
1009                 int const lines = sample_size(i).height;
1010                 for (int j = 0; j < lines; ++j) {
1011                         memcpy (p, q, _line_size[i]);
1012                         p += stride()[i];
1013                         q += other.stride()[i];
1014                 }
1015         }
1016 }
1017
1018 Image::Image (AVFrame* frame)
1019         : _size (frame->width, frame->height)
1020         , _pixel_format (static_cast<AVPixelFormat>(frame->format))
1021         , _aligned (true)
1022 {
1023         DCPOMATIC_ASSERT (_pixel_format != AV_PIX_FMT_NONE);
1024
1025         allocate ();
1026
1027         for (int i = 0; i < planes(); ++i) {
1028                 uint8_t* p = _data[i];
1029                 uint8_t* q = frame->data[i];
1030                 int const lines = sample_size(i).height;
1031                 for (int j = 0; j < lines; ++j) {
1032                         memcpy (p, q, _line_size[i]);
1033                         p += stride()[i];
1034                         /* AVFrame's linesize is what we call `stride' */
1035                         q += frame->linesize[i];
1036                 }
1037         }
1038 }
1039
1040 Image::Image (shared_ptr<const Image> other, bool aligned)
1041         : _size (other->_size)
1042         , _pixel_format (other->_pixel_format)
1043         , _aligned (aligned)
1044 {
1045         allocate ();
1046
1047         for (int i = 0; i < planes(); ++i) {
1048                 DCPOMATIC_ASSERT (line_size()[i] == other->line_size()[i]);
1049                 uint8_t* p = _data[i];
1050                 uint8_t* q = other->data()[i];
1051                 int const lines = sample_size(i).height;
1052                 for (int j = 0; j < lines; ++j) {
1053                         memcpy (p, q, line_size()[i]);
1054                         p += stride()[i];
1055                         q += other->stride()[i];
1056                 }
1057         }
1058 }
1059
1060 Image&
1061 Image::operator= (Image const & other)
1062 {
1063         if (this == &other) {
1064                 return *this;
1065         }
1066
1067         Image tmp (other);
1068         swap (tmp);
1069         return *this;
1070 }
1071
1072 void
1073 Image::swap (Image & other)
1074 {
1075         std::swap (_size, other._size);
1076         std::swap (_pixel_format, other._pixel_format);
1077
1078         for (int i = 0; i < 4; ++i) {
1079                 std::swap (_data[i], other._data[i]);
1080                 std::swap (_line_size[i], other._line_size[i]);
1081                 std::swap (_stride[i], other._stride[i]);
1082         }
1083
1084         std::swap (_aligned, other._aligned);
1085 }
1086
1087 Image::~Image ()
1088 {
1089         for (int i = 0; i < planes(); ++i) {
1090                 av_free (_data[i]);
1091         }
1092
1093         av_free (_data);
1094         av_free (_line_size);
1095         av_free (_stride);
1096 }
1097
1098 uint8_t * const *
1099 Image::data () const
1100 {
1101         return _data;
1102 }
1103
1104 int const *
1105 Image::line_size () const
1106 {
1107         return _line_size;
1108 }
1109
1110 int const *
1111 Image::stride () const
1112 {
1113         return _stride;
1114 }
1115
1116 dcp::Size
1117 Image::size () const
1118 {
1119         return _size;
1120 }
1121
1122 bool
1123 Image::aligned () const
1124 {
1125         return _aligned;
1126 }
1127
1128
1129 PositionImage
1130 merge (list<PositionImage> images)
1131 {
1132         if (images.empty ()) {
1133                 return {};
1134         }
1135
1136         if (images.size() == 1) {
1137                 return images.front ();
1138         }
1139
1140         dcpomatic::Rect<int> all (images.front().position, images.front().image->size().width, images.front().image->size().height);
1141         for (auto const& i: images) {
1142                 all.extend (dcpomatic::Rect<int>(i.position, i.image->size().width, i.image->size().height));
1143         }
1144
1145         auto merged = make_shared<Image>(images.front().image->pixel_format(), dcp::Size(all.width, all.height), true);
1146         merged->make_transparent ();
1147         for (auto const& i: images) {
1148                 merged->alpha_blend (i.image, i.position - all.position());
1149         }
1150
1151         return PositionImage (merged, all.position ());
1152 }
1153
1154
1155 bool
1156 operator== (Image const & a, Image const & b)
1157 {
1158         if (a.planes() != b.planes() || a.pixel_format() != b.pixel_format() || a.aligned() != b.aligned()) {
1159                 return false;
1160         }
1161
1162         for (int c = 0; c < a.planes(); ++c) {
1163                 if (a.sample_size(c).height != b.sample_size(c).height || a.line_size()[c] != b.line_size()[c] || a.stride()[c] != b.stride()[c]) {
1164                         return false;
1165                 }
1166
1167                 uint8_t* p = a.data()[c];
1168                 uint8_t* q = b.data()[c];
1169                 int const lines = a.sample_size(c).height;
1170                 for (int y = 0; y < lines; ++y) {
1171                         if (memcmp (p, q, a.line_size()[c]) != 0) {
1172                                 return false;
1173                         }
1174
1175                         p += a.stride()[c];
1176                         q += b.stride()[c];
1177                 }
1178         }
1179
1180         return true;
1181 }
1182
1183 /** Fade the image.
1184  *  @param f Amount to fade by; 0 is black, 1 is no fade.
1185  */
1186 void
1187 Image::fade (float f)
1188 {
1189         /* U/V black value for 8-bit colour */
1190         static int const eight_bit_uv =    (1 << 7) - 1;
1191         /* U/V black value for 10-bit colour */
1192         static uint16_t const ten_bit_uv = (1 << 9) - 1;
1193
1194         switch (_pixel_format) {
1195         case AV_PIX_FMT_YUV420P:
1196         {
1197                 /* Y */
1198                 uint8_t* p = data()[0];
1199                 int const lines = sample_size(0).height;
1200                 for (int y = 0; y < lines; ++y) {
1201                         uint8_t* q = p;
1202                         for (int x = 0; x < line_size()[0]; ++x) {
1203                                 *q = int(float(*q) * f);
1204                                 ++q;
1205                         }
1206                         p += stride()[0];
1207                 }
1208
1209                 /* U, V */
1210                 for (int c = 1; c < 3; ++c) {
1211                         uint8_t* p = data()[c];
1212                         int const lines = sample_size(c).height;
1213                         for (int y = 0; y < lines; ++y) {
1214                                 uint8_t* q = p;
1215                                 for (int x = 0; x < line_size()[c]; ++x) {
1216                                         *q = eight_bit_uv + int((int(*q) - eight_bit_uv) * f);
1217                                         ++q;
1218                                 }
1219                                 p += stride()[c];
1220                         }
1221                 }
1222
1223                 break;
1224         }
1225
1226         case AV_PIX_FMT_RGB24:
1227         {
1228                 /* 8-bit */
1229                 uint8_t* p = data()[0];
1230                 int const lines = sample_size(0).height;
1231                 for (int y = 0; y < lines; ++y) {
1232                         uint8_t* q = p;
1233                         for (int x = 0; x < line_size()[0]; ++x) {
1234                                 *q = int (float (*q) * f);
1235                                 ++q;
1236                         }
1237                         p += stride()[0];
1238                 }
1239                 break;
1240         }
1241
1242         case AV_PIX_FMT_XYZ12LE:
1243         case AV_PIX_FMT_RGB48LE:
1244                 /* 16-bit little-endian */
1245                 for (int c = 0; c < 3; ++c) {
1246                         int const stride_pixels = stride()[c] / 2;
1247                         int const line_size_pixels = line_size()[c] / 2;
1248                         uint16_t* p = reinterpret_cast<uint16_t*> (data()[c]);
1249                         int const lines = sample_size(c).height;
1250                         for (int y = 0; y < lines; ++y) {
1251                                 uint16_t* q = p;
1252                                 for (int x = 0; x < line_size_pixels; ++x) {
1253                                         *q = int (float (*q) * f);
1254                                         ++q;
1255                                 }
1256                                 p += stride_pixels;
1257                         }
1258                 }
1259                 break;
1260
1261         case AV_PIX_FMT_YUV422P10LE:
1262         {
1263                 /* Y */
1264                 {
1265                         int const stride_pixels = stride()[0] / 2;
1266                         int const line_size_pixels = line_size()[0] / 2;
1267                         uint16_t* p = reinterpret_cast<uint16_t*> (data()[0]);
1268                         int const lines = sample_size(0).height;
1269                         for (int y = 0; y < lines; ++y) {
1270                                 uint16_t* q = p;
1271                                 for (int x = 0; x < line_size_pixels; ++x) {
1272                                         *q = int(float(*q) * f);
1273                                         ++q;
1274                                 }
1275                                 p += stride_pixels;
1276                         }
1277                 }
1278
1279                 /* U, V */
1280                 for (int c = 1; c < 3; ++c) {
1281                         int const stride_pixels = stride()[c] / 2;
1282                         int const line_size_pixels = line_size()[c] / 2;
1283                         uint16_t* p = reinterpret_cast<uint16_t*> (data()[c]);
1284                         int const lines = sample_size(c).height;
1285                         for (int y = 0; y < lines; ++y) {
1286                                 uint16_t* q = p;
1287                                 for (int x = 0; x < line_size_pixels; ++x) {
1288                                         *q = ten_bit_uv + int((int(*q) - ten_bit_uv) * f);
1289                                         ++q;
1290                                 }
1291                                 p += stride_pixels;
1292                         }
1293                 }
1294                 break;
1295
1296         }
1297
1298         default:
1299                 throw PixelFormatError ("fade()", _pixel_format);
1300         }
1301 }
1302
1303 shared_ptr<const Image>
1304 Image::ensure_aligned (shared_ptr<const Image> image)
1305 {
1306         if (image->aligned()) {
1307                 return image;
1308         }
1309
1310         return make_shared<Image>(image, true);
1311 }
1312
1313 size_t
1314 Image::memory_used () const
1315 {
1316         size_t m = 0;
1317         for (int i = 0; i < planes(); ++i) {
1318                 m += _stride[i] * sample_size(i).height;
1319         }
1320         return m;
1321 }
1322
1323 class Memory
1324 {
1325 public:
1326         Memory ()
1327                 : data(0)
1328                 , size(0)
1329         {}
1330
1331         ~Memory ()
1332         {
1333                 free (data);
1334         }
1335
1336         uint8_t* data;
1337         size_t size;
1338 };
1339
1340 static void
1341 png_write_data (png_structp png_ptr, png_bytep data, png_size_t length)
1342 {
1343         auto mem = reinterpret_cast<Memory*>(png_get_io_ptr(png_ptr));
1344         size_t size = mem->size + length;
1345
1346         if (mem->data) {
1347                 mem->data = reinterpret_cast<uint8_t*>(realloc(mem->data, size));
1348         } else {
1349                 mem->data = reinterpret_cast<uint8_t*>(malloc(size));
1350         }
1351
1352         if (!mem->data) {
1353                 throw EncodeError (N_("could not allocate memory for PNG"));
1354         }
1355
1356         memcpy (mem->data + mem->size, data, length);
1357         mem->size += length;
1358 }
1359
1360 static void
1361 png_flush (png_structp)
1362 {
1363
1364 }
1365
1366 static void
1367 png_error_fn (png_structp png_ptr, char const * message)
1368 {
1369         reinterpret_cast<Image*>(png_get_error_ptr(png_ptr))->png_error (message);
1370 }
1371
1372 void
1373 Image::png_error (char const * message)
1374 {
1375         throw EncodeError (String::compose ("Error during PNG write: %1", message));
1376 }
1377
1378 dcp::ArrayData
1379 Image::as_png () const
1380 {
1381         DCPOMATIC_ASSERT (bytes_per_pixel(0) == 4);
1382         DCPOMATIC_ASSERT (planes() == 1);
1383         if (pixel_format() != AV_PIX_FMT_RGBA) {
1384                 return convert_pixel_format(dcp::YUVToRGB::REC709, AV_PIX_FMT_RGBA, true, false)->as_png();
1385         }
1386
1387         /* error handling? */
1388         png_structp png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING, reinterpret_cast<void*>(const_cast<Image*>(this)), png_error_fn, 0);
1389         if (!png_ptr) {
1390                 throw EncodeError (N_("could not create PNG write struct"));
1391         }
1392
1393         Memory state;
1394
1395         png_set_write_fn (png_ptr, &state, png_write_data, png_flush);
1396
1397         png_infop info_ptr = png_create_info_struct(png_ptr);
1398         if (!info_ptr) {
1399                 png_destroy_write_struct (&png_ptr, &info_ptr);
1400                 throw EncodeError (N_("could not create PNG info struct"));
1401         }
1402
1403         png_set_IHDR (png_ptr, info_ptr, size().width, size().height, 8, PNG_COLOR_TYPE_RGBA, PNG_INTERLACE_NONE, PNG_COMPRESSION_TYPE_DEFAULT, PNG_FILTER_TYPE_DEFAULT);
1404
1405         png_byte ** row_pointers = reinterpret_cast<png_byte **>(png_malloc(png_ptr, size().height * sizeof(png_byte *)));
1406         for (int i = 0; i < size().height; ++i) {
1407                 row_pointers[i] = (png_byte *) (data()[0] + i * stride()[0]);
1408         }
1409
1410         png_write_info (png_ptr, info_ptr);
1411         png_write_image (png_ptr, row_pointers);
1412         png_write_end (png_ptr, info_ptr);
1413
1414         png_destroy_write_struct (&png_ptr, &info_ptr);
1415         png_free (png_ptr, row_pointers);
1416
1417         return dcp::ArrayData (state.data, state.size);
1418 }
1419
1420
1421 void
1422 Image::video_range_to_full_range ()
1423 {
1424         switch (_pixel_format) {
1425         case AV_PIX_FMT_RGB24:
1426         {
1427                 float const factor = 256.0 / 219.0;
1428                 uint8_t* p = data()[0];
1429                 int const lines = sample_size(0).height;
1430                 for (int y = 0; y < lines; ++y) {
1431                         uint8_t* q = p;
1432                         for (int x = 0; x < line_size()[0]; ++x) {
1433                                 *q = int((*q - 16) * factor);
1434                                 ++q;
1435                         }
1436                         p += stride()[0];
1437                 }
1438                 break;
1439         }
1440         case AV_PIX_FMT_GBRP12LE:
1441         {
1442                 float const factor = 4096.0 / 3504.0;
1443                 for (int c = 0; c < 3; ++c) {
1444                         uint16_t* p = reinterpret_cast<uint16_t*>(data()[c]);
1445                         int const lines = sample_size(c).height;
1446                         for (int y = 0; y < lines; ++y) {
1447                                 uint16_t* q = p;
1448                                 int const line_size_pixels = line_size()[c] / 2;
1449                                 for (int x = 0; x < line_size_pixels; ++x) {
1450                                         *q = int((*q - 256) * factor);
1451                                         ++q;
1452                                 }
1453                         }
1454                 }
1455                 break;
1456         }
1457         default:
1458                 throw PixelFormatError ("video_range_to_full_range()", _pixel_format);
1459         }
1460 }
1461