GNU Octave  4.0.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
__magick_read__.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 2013-2015 CarnĂ« Draug
4 Copyright (C) 2002-2015 Andy Adler
5 Copyright (C) 2008 Thomas L. Scofield
6 Copyright (C) 2010 David Grundberg
7 
8 This file is part of Octave.
9 
10 Octave is free software; you can redistribute it and/or modify it
11 under the terms of the GNU General Public License as published by the
12 Free Software Foundation; either version 3 of the License, or (at your
13 option) any later version.
14 
15 Octave is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
17 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
18 for more details.
19 
20 You should have received a copy of the GNU General Public License
21 along with Octave; see the file COPYING. If not, see
22 <http://www.gnu.org/licenses/>.
23 
24 */
25 
26 #ifdef HAVE_CONFIG_H
27 #include <config.h>
28 #endif
29 
30 #include "file-stat.h"
31 #include "oct-env.h"
32 #include "oct-time.h"
33 
34 #include "defun-dld.h"
35 #include "error.h"
36 #include "ov-struct.h"
37 
38 #include "gripes.h"
39 
40 #ifdef HAVE_MAGICK
41 
42 #include <Magick++.h>
43 #include <clocale>
44 
45 // In theory, it should be enough to check the class:
46 // Magick::ClassType
47 // PseudoClass:
48 // Image is composed of pixels which specify an index in a color palette.
49 // DirectClass:
50 // Image is composed of pixels which represent literal color values.
51 //
52 // GraphicsMagick does not really distinguishes between indexed and
53 // normal images. After reading a file, it decides itself the optimal
54 // way to store the image in memory, independently of the how the
55 // image was stored in the file. That's what ClassType returns. While
56 // it seems to match the original file most of the times, this is
57 // not necessarily true all the times. See
58 // https://sourceforge.net/mailarchive/message.php?msg_id=31180507
59 // In addition to the ClassType, there is also ImageType which has a
60 // type for indexed images (PaletteType and PaletteMatteType). However,
61 // they also don't represent the original image. Not only does DirectClass
62 // can have a PaletteType, but also does a PseudoClass have non Palette
63 // types.
64 //
65 // We can't do better without having format specific code which is
66 // what we are trying to avoid by using a library such as GM. We at
67 // least create workarounds for the most common problems.
68 //
69 // 1) A grayscale jpeg image can report being indexed even though the
70 // JPEG format has no support for indexed images. We can at least
71 // fix this one.
72 // 2) A PNG file is only an indexed image if color type orig is 3 (value comes
73 // from libpng)
74 static bool
75 is_indexed (const Magick::Image& img)
76 {
77  bool indexed = (img.classType () == Magick::PseudoClass);
78  // Our problem until now is non-indexed images, being represented as indexed
79  // by GM. The following attempts educated guesses to undo this optimization.
80  if (indexed)
81  {
82  const std::string fmt = img.magick ();
83  if (fmt == "JPEG")
84  // The JPEG format does not support indexed images, but GM sometimes
85  // reports grayscale JPEG as indexed. Always false for JPEG.
86  indexed = false;
87  else if (fmt == "PNG")
88  {
89  // Newer versions of GM (at least does not happens with 1.3.16) will
90  // store values from the underlying library as image attributes. In
91  // the case of PNG files, this is libpng where an indexed image will
92  // always have a value of 3 for "color-type-orig". This property
93  // always has a value in libpng so if we get nothing, we assume this
94  // GM version does not store them and we have to go with whatever
95  // GM PseudoClass says.
96  const std::string color_type =
97  const_cast<Magick::Image&> (img).attribute ("PNG:IHDR.color-type-orig");
98  if (! color_type.empty() && color_type != "3")
99  indexed = false;
100  }
101  }
102  return indexed;
103 }
104 
105 // The depth from depth() is not always correct for us but seems to be the
106 // best value we can get. For example, a grayscale png image with 1 bit
107 // per channel should return a depth of 1 but instead we get 8.
108 // We could check channelDepth() but then, which channel has the data
109 // is not straightforward. So we'd have to check all
110 // the channels and select the highest value. But then, I also
111 // have a 16bit TIFF whose depth returns 16 (correct), but all of the
112 // channels gives 8 (wrong). No idea why, maybe a bug in GM?
113 // Anyway, using depth() seems that only causes problems for binary
114 // images, and the problem with channelDepth() is not making set them
115 // all to 1. So we will guess that if all channels have depth of 1,
116 // then we must have a binary image.
117 // Note that we can't use AllChannels it doesn't work for this.
118 // Instead of checking all of the individual channels, we check one
119 // from RGB, CMYK, grayscale, and transparency.
120 static octave_idx_type
121 get_depth (Magick::Image& img)
122 {
123  octave_idx_type depth = img.depth ();
124  if (depth == 8
125  && img.channelDepth (Magick::RedChannel) == 1
126  && img.channelDepth (Magick::CyanChannel) == 1
127  && img.channelDepth (Magick::OpacityChannel) == 1
128  && img.channelDepth (Magick::GrayChannel) == 1)
129  depth = 1;
130 
131  return depth;
132 }
133 
134 // We need this in case one of the sides of the image being read has
135 // width 1. In those cases, the type will come as scalar instead of range
136 // since that's the behaviour of the colon operator (1:1:1 will be a scalar,
137 // not a range).
138 static Range
140 {
141  Range output;
142  if (region.is_range ())
143  output = region.range_value ();
144  else if (region.is_scalar_type ())
145  {
146  double value = region.scalar_value ();
147  output = Range (value, value);
148  }
149  else
150  error ("__magick_read__: unknow datatype for Region option");
151 
152  return output;
153 }
154 
155 static std::map<std::string, octave_idx_type>
157 {
158  std::map<std::string, octave_idx_type> region;
159  const Cell pixel_region = options.getfield ("region").cell_value ();
160 
161  // Subtract 1 to account for 0 indexing.
162  const Range rows = get_region_range (pixel_region (0));
163  const Range cols = get_region_range (pixel_region (1));
164  region["row_start"] = rows.base () -1;
165  region["col_start"] = cols.base () -1;
166  region["row_end"] = rows.max () -1;
167  region["col_end"] = cols.max () -1;
168 
169  // Length of the area to load into the Image Pixel Cache. We use max and
170  // min to account for cases where last element of range is the range limit.
171  region["row_cache"] = region["row_end"] - region["row_start"] +1;
172  region["col_cache"] = region["col_end"] - region["col_start"] +1;
173 
174  // How much we have to shift in the memory when doing the loops.
175  region["row_shift"] = region["col_cache"] * rows.inc ();
176  region["col_shift"] = region["col_cache"] *
177  (region["row_cache"] + rows.inc () -1) - cols.inc ();
178 
179  // The actual height and width of the output image
180  region["row_out"] = rows.nelem ();
181  region["col_out"] = cols.nelem ();
182 
183  return region;
184 }
185 
186 static octave_value_list
187 read_maps (Magick::Image& img)
188 {
189  // can't call colorMapSize on const Magick::Image
190  const octave_idx_type mapsize = img.colorMapSize ();
191  Matrix cmap = Matrix (mapsize, 3); // colormap
192  ColumnVector amap = ColumnVector (mapsize); // alpha map
193  for (octave_idx_type i = 0; i < mapsize; i++)
194  {
195  const Magick::ColorRGB c = img.colorMap (i);
196  cmap(i,0) = c.red ();
197  cmap(i,1) = c.green ();
198  cmap(i,2) = c.blue ();
199  amap(i) = c.alpha ();
200  }
201  octave_value_list maps;
202  maps(0) = cmap;
203  maps(1) = amap;
204  return maps;
205 }
206 
207 template <class T>
208 static octave_value_list
209 read_indexed_images (const std::vector<Magick::Image>& imvec,
210  const Array<octave_idx_type>& frameidx,
211  const octave_idx_type& nargout,
212  const octave_scalar_map& options)
213 {
214  typedef typename T::element_type P;
215 
216  octave_value_list retval (3, Matrix ());
217 
218  std::map<std::string, octave_idx_type> region = calculate_region (options);
219  const octave_idx_type nFrames = frameidx.length ();
220  const octave_idx_type nRows = region["row_out"];
221  const octave_idx_type nCols = region["col_out"];
222 
223  // imvec has all of the pages of a file, even the ones we are not
224  // interested in. We will use the first image that we will be actually
225  // reading to get information about the image.
226  const octave_idx_type def_elem = frameidx(0);
227 
228  T img = T (dim_vector (nRows, nCols, 1, nFrames));
229  P* img_fvec = img.fortran_vec ();
230 
231  const octave_idx_type row_start = region["row_start"];
232  const octave_idx_type col_start = region["col_start"];
233  const octave_idx_type row_shift = region["row_shift"];
234  const octave_idx_type col_shift = region["col_shift"];
235  const octave_idx_type row_cache = region["row_cache"];
236  const octave_idx_type col_cache = region["col_cache"];
237 
238  // When reading PixelPackets from the Image Pixel Cache, they come in
239  // row major order. So we keep moving back and forth there so we can
240  // write the image in column major order.
241  octave_idx_type idx = 0;
242  for (octave_idx_type frame = 0; frame < nFrames; frame++)
243  {
244  OCTAVE_QUIT;
245  imvec[frameidx(frame)].getConstPixels (col_start, row_start,
246  col_cache, row_cache);
247 
248  const Magick::IndexPacket *pix
249  = imvec[frameidx(frame)].getConstIndexes ();
250 
251  for (octave_idx_type col = 0; col < nCols; col++)
252  {
253  for (octave_idx_type row = 0; row < nRows; row++)
254  {
255  img_fvec[idx++] = static_cast<P> (*pix);
256  pix += row_shift;
257  }
258  pix -= col_shift;
259  }
260  }
261  retval(0) = octave_value (img);
262 
263  // Only bother reading the colormap if it was requested as output.
264  if (nargout > 1)
265  {
266  // In theory, it should be possible for each frame of an image to
267  // have different colormaps but for Matlab compatibility, we only
268  // return the colormap of the first frame. To obtain the colormaps
269  // of different frames, one needs can either use imfinfo or a for
270  // loop around imread.
271  const octave_value_list maps =
272  read_maps (const_cast<Magick::Image&> (imvec[frameidx(def_elem)]));
273 
274  retval(1) = maps(0);
275 
276  // only interpret alpha channel if it exists and was requested as output
277  if (imvec[def_elem].matte () && nargout >= 3)
278  {
279  const Matrix amap = maps(1).matrix_value ();
280  const double* amap_fvec = amap.fortran_vec ();
281 
282  NDArray alpha (dim_vector (nRows, nCols, 1, nFrames));
283  double* alpha_fvec = alpha.fortran_vec ();
284 
285  // GraphicsMagick stores the alpha values inverted, i.e.,
286  // 1 for transparent and 0 for opaque so we fix that here.
287  const octave_idx_type nPixels = alpha.numel ();
288  for (octave_idx_type pix = 0; pix < nPixels; pix++)
289  alpha_fvec[pix] = 1 - amap_fvec[static_cast<int> (img_fvec[3])];
290 
291  retval(2) = alpha;
292  }
293  }
294 
295  return retval;
296 }
297 
298 // This function is highly repetitive, a bunch of for loops that are
299 // very similar to account for different image types. They are different
300 // enough that trying to reduce the copy and paste would decrease its
301 // readability too much.
302 template <class T>
304 read_images (std::vector<Magick::Image>& imvec,
305  const Array<octave_idx_type>& frameidx,
306  const octave_idx_type& nargout,
307  const octave_scalar_map& options)
308 {
309  typedef typename T::element_type P;
310 
311  octave_value_list retval (3, Matrix ());
312 
313  std::map<std::string, octave_idx_type> region = calculate_region (options);
314  const octave_idx_type nFrames = frameidx.length ();
315  const octave_idx_type nRows = region["row_out"];
316  const octave_idx_type nCols = region["col_out"];
317  T img;
318 
319  // imvec has all of the pages of a file, even the ones we are not
320  // interested in. We will use the first image that we will be actually
321  // reading to get information about the image.
322  const octave_idx_type def_elem = frameidx(0);
323 
324  const octave_idx_type row_start = region["row_start"];
325  const octave_idx_type col_start = region["col_start"];
326  const octave_idx_type row_shift = region["row_shift"];
327  const octave_idx_type col_shift = region["col_shift"];
328  const octave_idx_type row_cache = region["row_cache"];
329  const octave_idx_type col_cache = region["col_cache"];
330 
331  // GraphicsMagick (GM) keeps the image values in memory using whatever
332  // QuantumDepth it was built with independently of the original image
333  // bitdepth. Basically this means that if GM was built with quantum 16
334  // all values are scaled in the uint16 range. If the original image
335  // had an 8 bit depth, we need to rescale it for that range.
336  // However, if the image had a bitdepth of 32, then we will be returning
337  // a floating point image. In this case, the values need to be rescaled
338  // for the range [0 1] (this is what Matlab has documented on the page
339  // about image types but in some cases seems to be doing something else.
340  // See bug #39249).
341  // Finally, we must do the division ourselves (set a divisor) instead of
342  // using quantumOperator for the cases where we will be returning floating
343  // point and want things in the range [0 1]. This is the same reason why
344  // the divisor is of type double.
345  // uint64_t is used in expression because default 32-bit value overflows
346  // when depth() is 32.
347  // TODO in the next release of GraphicsMagick, MaxRGB should be replaced
348  // with QuantumRange since MaxRGB is already deprecated in ImageMagick.
349  double divisor;
350  if (imvec[def_elem].depth () == 32)
352  else
353  divisor = MaxRGB / ((uint64_t (1) << imvec[def_elem].depth ()) - 1);
354 
355  // FIXME: this workaround should probably be fixed in GM by creating a
356  // new ImageType BilevelMatteType
357  // Despite what GM documentation claims, opacity is not only on the types
358  // with Matte on the name. It is possible that an image is completely
359  // black (1 color), and have a second channel set for transparency (2nd
360  // color). Its type will be bilevel since there is no BilevelMatte. The
361  // only way to check for this seems to be by checking matte ().
362  Magick::ImageType type = imvec[def_elem].type ();
363  if (type == Magick::BilevelType && imvec[def_elem].matte ())
364  type = Magick::GrayscaleMatteType;
365 
366  // FIXME: ImageType is the type being used to represent the image in memory
367  // by GM. The real type may be different (see among others bug #36820). For
368  // example, a png file where all channels are equal may report being
369  // grayscale or even bilevel. But we must always return the real image in
370  // file. In some cases, the original image attributes are stored in the
371  // attributes but this is undocumented. This should be fixed in GM so that
372  // a method such as original_type returns an actual Magick::ImageType
373  if (imvec[0].magick () == "PNG")
374  {
375  // These values come from libpng, not GM:
376  // Grayscale = 0
377  // Palette = 2 + 1
378  // RGB = 2
379  // RGB + Alpha = 2 + 4
380  // Grayscale + Alpha = 4
381  // We won't bother with case 3 (palette) since those should be
382  // read by the function to read indexed images
383  const std::string type_str
384  = imvec[0].attribute ("PNG:IHDR.color-type-orig");
385 
386  if (type_str == "0")
387  type = Magick::GrayscaleType;
388  else if (type_str == "2")
389  type = Magick::TrueColorType;
390  else if (type_str == "6")
391  type = Magick::TrueColorMatteType;
392  else if (type_str == "4")
393  type = Magick::GrayscaleMatteType;
394  // Color types 0, 2, and 3 can also have alpha channel, conveyed
395  // via the "tRNS" chunk. For 0 and 2, it's limited to GIF-style
396  // binary transparency, while 3 can have any level of alpha per
397  // palette entry. We thus must check matte() to see if the image
398  // really doesn't have an alpha channel.
399  if (imvec[0].matte ())
400  {
401  if (type == Magick::GrayscaleType)
402  type = Magick::GrayscaleMatteType;
403  else if (type == Magick::TrueColorType)
404  type = Magick::TrueColorMatteType;
405  }
406  }
407 
408  // If the alpha channel was not requested, treat images as if
409  // it doesn't exist.
410  if (nargout < 3)
411  {
412  switch (type)
413  {
414  case Magick::GrayscaleMatteType:
415  type = Magick::GrayscaleType;
416  break;
417 
418  case Magick::PaletteMatteType:
419  type = Magick::PaletteType;
420  break;
421 
422  case Magick::TrueColorMatteType:
423  type = Magick::TrueColorType;
424  break;
425 
426  case Magick::ColorSeparationMatteType:
427  type = Magick::ColorSeparationType;
428  break;
429 
430  default:
431  // Do nothing other than silencing warnings about enumeration
432  // values not being handled in switch.
433  ;
434  }
435  }
436 
437  const octave_idx_type colour_stride = nRows * nCols;
438  switch (type)
439  {
440  case Magick::BilevelType: // Monochrome bi-level image
441  case Magick::GrayscaleType: // Grayscale image
442  {
443  img = T (dim_vector (nRows, nCols, 1, nFrames));
444  P *img_fvec = img.fortran_vec ();
445 
446  octave_idx_type idx = 0;
447  for (octave_idx_type frame = 0; frame < nFrames; frame++)
448  {
449  OCTAVE_QUIT;
450  const Magick::PixelPacket *pix
451  = imvec[frameidx(frame)].getConstPixels (col_start, row_start,
452  col_cache, row_cache);
453 
454  for (octave_idx_type col = 0; col < nCols; col++)
455  {
456  for (octave_idx_type row = 0; row < nRows; row++)
457  {
458  img_fvec[idx++] = pix->red / divisor;
459  pix += row_shift;
460  }
461  pix -= col_shift;
462  }
463  }
464  break;
465  }
466 
467  case Magick::GrayscaleMatteType: // Grayscale image with opacity
468  {
469  img = T (dim_vector (nRows, nCols, 1, nFrames));
470  T alpha (dim_vector (nRows, nCols, 1, nFrames));
471  P *img_fvec = img.fortran_vec ();
472  P *a_fvec = alpha.fortran_vec ();
473 
474  octave_idx_type idx = 0;
475  for (octave_idx_type frame = 0; frame < nFrames; frame++)
476  {
477  OCTAVE_QUIT;
478  const Magick::PixelPacket *pix
479  = imvec[frameidx(frame)].getConstPixels (col_start, row_start,
480  col_cache, row_cache);
481 
482  for (octave_idx_type col = 0; col < nCols; col++)
483  {
484  for (octave_idx_type row = 0; row < nRows; row++)
485  {
486  img_fvec[idx] = pix->red / divisor;
487  a_fvec[idx] = (MaxRGB - pix->opacity) / divisor;
488  pix += row_shift;
489  idx++;
490  }
491  pix -= col_shift;
492  }
493  }
494  retval(2) = alpha;
495  break;
496  }
497 
498  case Magick::PaletteType: // Indexed color (palette) image
499  case Magick::TrueColorType: // Truecolor image
500  {
501  img = T (dim_vector (nRows, nCols, 3, nFrames));
502  P *img_fvec = img.fortran_vec ();
503 
504  const octave_idx_type frame_stride = colour_stride * 3;
505  for (octave_idx_type frame = 0; frame < nFrames; frame++)
506  {
507  OCTAVE_QUIT;
508  const Magick::PixelPacket *pix
509  = imvec[frameidx(frame)].getConstPixels (col_start, row_start,
510  col_cache, row_cache);
511 
512  octave_idx_type idx = 0;
513  P *rbuf = img_fvec;
514  P *gbuf = img_fvec + colour_stride;
515  P *bbuf = img_fvec + colour_stride * 2;
516 
517  for (octave_idx_type col = 0; col < nCols; col++)
518  {
519  for (octave_idx_type row = 0; row < nRows; row++)
520  {
521  rbuf[idx] = pix->red / divisor;
522  gbuf[idx] = pix->green / divisor;
523  bbuf[idx] = pix->blue / divisor;
524  pix += row_shift;
525  idx++;
526  }
527  pix -= col_shift;
528  }
529  img_fvec += frame_stride;
530  }
531  break;
532  }
533 
534  case Magick::PaletteMatteType: // Indexed color image with opacity
535  case Magick::TrueColorMatteType: // Truecolor image with opacity
536  {
537  img = T (dim_vector (nRows, nCols, 3, nFrames));
538  T alpha (dim_vector (nRows, nCols, 1, nFrames));
539  P *img_fvec = img.fortran_vec ();
540  P *a_fvec = alpha.fortran_vec ();
541 
542  const octave_idx_type frame_stride = colour_stride * 3;
543 
544  // Unlike the index for the other channels, this one won't need
545  // to be reset on each frame since it's a separate matrix.
546  octave_idx_type a_idx = 0;
547  for (octave_idx_type frame = 0; frame < nFrames; frame++)
548  {
549  OCTAVE_QUIT;
550  const Magick::PixelPacket *pix
551  = imvec[frameidx(frame)].getConstPixels (col_start, row_start,
552  col_cache, row_cache);
553 
554  octave_idx_type idx = 0;
555  P *rbuf = img_fvec;
556  P *gbuf = img_fvec + colour_stride;
557  P *bbuf = img_fvec + colour_stride * 2;
558 
559  for (octave_idx_type col = 0; col < nCols; col++)
560  {
561  for (octave_idx_type row = 0; row < nRows; row++)
562  {
563  rbuf[idx] = pix->red / divisor;
564  gbuf[idx] = pix->green / divisor;
565  bbuf[idx] = pix->blue / divisor;
566  a_fvec[a_idx++] = (MaxRGB - pix->opacity) / divisor;
567  pix += row_shift;
568  idx++;
569  }
570  pix -= col_shift;
571  }
572  img_fvec += frame_stride;
573  }
574  retval(2) = alpha;
575  break;
576  }
577 
578  case Magick::ColorSeparationType: // Cyan/Magenta/Yellow/Black (CMYK) image
579  {
580  img = T (dim_vector (nRows, nCols, 4, nFrames));
581  P *img_fvec = img.fortran_vec ();
582 
583  const octave_idx_type frame_stride = colour_stride * 4;
584  for (octave_idx_type frame = 0; frame < nFrames; frame++)
585  {
586  OCTAVE_QUIT;
587  const Magick::PixelPacket *pix
588  = imvec[frameidx(frame)].getConstPixels (col_start, row_start,
589  col_cache, row_cache);
590 
591  octave_idx_type idx = 0;
592  P *cbuf = img_fvec;
593  P *mbuf = img_fvec + colour_stride;
594  P *ybuf = img_fvec + colour_stride * 2;
595  P *kbuf = img_fvec + colour_stride * 3;
596 
597  for (octave_idx_type col = 0; col < nCols; col++)
598  {
599  for (octave_idx_type row = 0; row < nRows; row++)
600  {
601  cbuf[idx] = pix->red / divisor;
602  mbuf[idx] = pix->green / divisor;
603  ybuf[idx] = pix->blue / divisor;
604  kbuf[idx] = pix->opacity / divisor;
605  pix += row_shift;
606  idx++;
607  }
608  pix -= col_shift;
609  }
610  img_fvec += frame_stride;
611  }
612  break;
613  }
614 
615  // Cyan, magenta, yellow, and black with alpha (opacity) channel
616  case Magick::ColorSeparationMatteType:
617  {
618  img = T (dim_vector (nRows, nCols, 4, nFrames));
619  T alpha (dim_vector (nRows, nCols, 1, nFrames));
620  P *img_fvec = img.fortran_vec ();
621  P *a_fvec = alpha.fortran_vec ();
622 
623  const octave_idx_type frame_stride = colour_stride * 4;
624 
625  // Unlike the index for the other channels, this one won't need
626  // to be reset on each frame since it's a separate matrix.
627  octave_idx_type a_idx = 0;
628  for (octave_idx_type frame = 0; frame < nFrames; frame++)
629  {
630  OCTAVE_QUIT;
631  const Magick::PixelPacket *pix
632  = imvec[frameidx(frame)].getConstPixels (col_start, row_start,
633  col_cache, row_cache);
634  // Note that for CMYKColorspace + matte (CMYKA), the opacity is
635  // stored in the assocated IndexPacket.
636  const Magick::IndexPacket *apix
637  = imvec[frameidx(frame)].getConstIndexes ();
638 
639  octave_idx_type idx = 0;
640  P *cbuf = img_fvec;
641  P *mbuf = img_fvec + colour_stride;
642  P *ybuf = img_fvec + colour_stride * 2;
643  P *kbuf = img_fvec + colour_stride * 3;
644 
645  for (octave_idx_type col = 0; col < nCols; col++)
646  {
647  for (octave_idx_type row = 0; row < nRows; row++)
648  {
649  cbuf[idx] = pix->red / divisor;
650  mbuf[idx] = pix->green / divisor;
651  ybuf[idx] = pix->blue / divisor;
652  kbuf[idx] = pix->opacity / divisor;
653  a_fvec[a_idx++] = (MaxRGB - *apix) / divisor;
654  pix += row_shift;
655  idx++;
656  }
657  pix -= col_shift;
658  }
659  img_fvec += frame_stride;
660  }
661  retval(2) = alpha;
662  break;
663  }
664 
665  default:
666  error ("__magick_read__: unknown Magick++ image type");
667  return retval;
668  }
669 
670  retval(0) = img;
671  return retval;
672 }
673 
674 // Read a file into vector of image objects.
675 void static
676 read_file (const std::string& filename, std::vector<Magick::Image>& imvec)
677 {
678  try
679  {
680  Magick::readImages (&imvec, filename);
681  }
682  catch (Magick::Warning& w)
683  {
684  warning ("Magick++ warning: %s", w.what ());
685  }
686  catch (Magick::Exception& e)
687  {
688  error ("Magick++ exception: %s", e.what ());
689  error_state = 1;
690  }
691 }
692 
693 static void
695 {
696  static bool initialized = false;
697 
698  if (! initialized)
699  {
700  // Save locale as GraphicsMagick might change this (fixed in
701  // GraphicsMagick since version 1.3.13 released on December 24, 2011)
702  const char *static_locale = setlocale (LC_ALL, NULL);
703  const std::string locale (static_locale);
704 
705  const std::string program_name
707  Magick::InitializeMagick (program_name.c_str ());
708 
709  // Restore locale from before GraphicsMagick initialisation
710  setlocale (LC_ALL, locale.c_str ());
711 
712  if (QuantumDepth < 32)
713  warning_with_id ("Octave:GraphicsMagic-Quantum-Depth",
714  "your version of %s limits images to %d bits per pixel",
715  MagickPackageName, QuantumDepth);
716 
717  initialized = true;
718  }
719 }
720 #endif
721 
722 DEFUN_DLD (__magick_read__, args, nargout,
723  "-*- texinfo -*-\n\
724 @deftypefn {Loadable Function} {[@var{img}, @var{map}, @var{alpha}] =} __magick_read__ (@var{fname}, @var{options})\n\
725 Read image with GraphicsMagick or ImageMagick.\n\
726 \n\
727 This is a private internal function not intended for direct use.\n\
728 Use @code{imread} instead.\n\
729 \n\
730 @seealso{imfinfo, imformats, imread, imwrite}\n\
731 @end deftypefn")
732 {
733  octave_value_list output;
734 
735 #ifndef HAVE_MAGICK
736  gripe_disabled_feature ("imread", "Image IO");
737 #else
738 
740 
741  if (args.length () != 2 || ! args(0).is_string ())
742  {
743  print_usage ();
744  return output;
745  }
746 
747  const octave_scalar_map options = args(1).scalar_map_value ();
748  if (error_state)
749  {
750  error ("__magick_read__: OPTIONS must be a struct");
751  return output;
752  }
753 
754  std::vector<Magick::Image> imvec;
755  read_file (args(0).string_value (), imvec);
756  if (error_state)
757  return output;
758 
759  // Prepare an Array with the indexes for the requested frames.
760  const octave_idx_type nFrames = imvec.size ();
761  Array<octave_idx_type> frameidx;
762  const octave_value indexes = options.getfield ("index");
763  if (indexes.is_string () && indexes.string_value () == "all")
764  {
765  frameidx.resize (dim_vector (1, nFrames));
766  for (octave_idx_type i = 0; i < nFrames; i++)
767  frameidx(i) = i;
768  }
769  else
770  {
771  frameidx = indexes.int_vector_value ();
772  if (error_state)
773  {
774  error ("__magick_read__: invalid value for Index/Frame");
775  return output;
776  }
777  // Fix indexes from base 1 to base 0, and at the same time, make
778  // sure none of the indexes is outside the range of image number.
779  const octave_idx_type n = frameidx.nelem ();
780  for (octave_idx_type i = 0; i < n; i++)
781  {
782  frameidx(i)--;
783  if (frameidx(i) < 0 || frameidx(i) > nFrames - 1)
784  {
785  // We do this check inside the loop because frameidx does not
786  // need to be ordered (this is a feature and even allows for
787  // some frames to be read multiple times).
788  error ("imread: index/frames specified are outside the number of images");
789  return output;
790  }
791  }
792  }
793 
794  // Check that all frames have the same size. We don't do this at the same
795  // time we decode the image because that's done in many different places,
796  // to cover the different types of images which would lead to a lot of
797  // copy and paste.
798  {
799  const unsigned int nRows = imvec[frameidx(0)].rows ();
800  const unsigned int nCols = imvec[frameidx(0)].columns ();
801  const octave_idx_type n = frameidx.nelem ();
802  for (octave_idx_type frame = 0; frame < n; frame++)
803  {
804  if (nRows != imvec[frameidx(frame)].rows ()
805  || nCols != imvec[frameidx(frame)].columns ())
806  {
807  error ("imread: all frames must have the same size but frame %i is different",
808  frameidx(frame) +1);
809  return output;
810  }
811  }
812  }
813 
814  const octave_idx_type depth = get_depth (imvec[frameidx(0)]);
815  if (is_indexed (imvec[frameidx(0)]))
816  {
817  if (depth <= 1)
818  output = read_indexed_images<boolNDArray> (imvec, frameidx,
819  nargout, options);
820  else if (depth <= 8)
821  output = read_indexed_images<uint8NDArray> (imvec, frameidx,
822  nargout, options);
823  else if (depth <= 16)
824  output = read_indexed_images<uint16NDArray> (imvec, frameidx,
825  nargout, options);
826  else
827  {
828  error ("imread: indexed images with depths greater than 16-bit are not supported");
829  return output;
830  }
831  }
832 
833  else
834  {
835  if (depth <= 1)
836  output = read_images<boolNDArray> (imvec, frameidx, nargout, options);
837  else if (depth <= 8)
838  output = read_images<uint8NDArray> (imvec, frameidx, nargout, options);
839  else if (depth <= 16)
840  output = read_images<uint16NDArray> (imvec, frameidx, nargout, options);
841  else if (depth <= 32)
842  output = read_images<FloatNDArray> (imvec, frameidx, nargout, options);
843  else
844  {
845  error ("imread: reading of images with %i-bit depth is not supported",
846  depth);
847  }
848  }
849 
850 #endif
851  return output;
852 }
853 
854 /*
855 ## No test needed for internal helper function.
856 %!assert (1)
857 */
858 
859 #ifdef HAVE_MAGICK
860 
861 template <class T>
862 static uint32NDArray
863 img_float2uint (const T& img)
864 {
865  typedef typename T::element_type P;
866  uint32NDArray out (img.dims ());
867 
868  octave_uint32* out_fvec = out.fortran_vec ();
869  const P* img_fvec = img.fortran_vec ();
870 
872  const octave_idx_type numel = img.numel ();
873  for (octave_idx_type idx = 0; idx < numel; idx++)
874  out_fvec[idx] = img_fvec[idx] * max;
875 
876  return out;
877 }
878 
879 // Gets the bitdepth to be used for an Octave class, i.e, returns 8 for
880 // uint8, 16 for uint16, and 32 for uint32
881 template <class T>
882 static octave_idx_type
884 {
885  typedef typename T::element_type P;
886  const octave_idx_type bitdepth =
887  sizeof (P) * std::numeric_limits<unsigned char>::digits;
888  return bitdepth;
889 }
890 
891 static Magick::Image
893  const octave_idx_type& bitdepth,
894  const Magick::ImageType& type,
895  const Magick::ClassType& klass)
896 {
897  Magick::Image img (Magick::Geometry (nCols, nRows), "black");
898  // Ensure that there are no other references to this image.
899  img.modifyImage ();
900 
901  img.classType (klass);
902  img.type (type);
903  // FIXME: for some reason, setting bitdepth doesn't seem to work for
904  // indexed images.
905  img.depth (bitdepth);
906  switch (type)
907  {
908  case Magick::GrayscaleMatteType:
909  case Magick::TrueColorMatteType:
910  case Magick::ColorSeparationMatteType:
911  case Magick::PaletteMatteType:
912  img.matte (true);
913  break;
914 
915  default:
916  img.matte (false);
917  }
918 
919  return img;
920 }
921 
922 template <class T>
923 static void
924 encode_indexed_images (std::vector<Magick::Image>& imvec,
925  const T& img,
926  const Matrix& cmap)
927 {
928  typedef typename T::element_type P;
929  const octave_idx_type nFrames = img.ndims () < 4 ? 1 : img.dims ()(3);
930  const octave_idx_type nRows = img.rows ();
931  const octave_idx_type nCols = img.columns ();
932  const octave_idx_type cmap_size = cmap.rows ();
933  const octave_idx_type bitdepth = bitdepth_from_class<T> ();
934 
935  // There is no colormap object, we need to build a new one for each frame,
936  // even if it's always the same. We can least get a vector for the Colors.
937  std::vector<Magick::ColorRGB> colormap;
938  {
939  const double* cmap_fvec = cmap.fortran_vec ();
940  const octave_idx_type G_offset = cmap_size;
941  const octave_idx_type B_offset = cmap_size * 2;
942  for (octave_idx_type map_idx = 0; map_idx < cmap_size; map_idx++)
943  colormap.push_back (Magick::ColorRGB (cmap_fvec[map_idx],
944  cmap_fvec[map_idx + G_offset],
945  cmap_fvec[map_idx + B_offset]));
946  }
947 
948  for (octave_idx_type frame = 0; frame < nFrames; frame++)
949  {
950  OCTAVE_QUIT;
951  Magick::Image m_img = init_enconde_image (nCols, nRows, bitdepth,
952  Magick::PaletteType,
953  Magick::PseudoClass);
954 
955  // Insert colormap.
956  m_img.colorMapSize (cmap_size);
957  for (octave_idx_type map_idx = 0; map_idx < cmap_size; map_idx++)
958  m_img.colorMap (map_idx, colormap[map_idx]);
959 
960  // Why are we also setting the pixel values instead of only the
961  // index values? We don't know if a file format supports indexed
962  // images. If we only set the indexes and then try to save the
963  // image as JPEG for example, the indexed values get discarded,
964  // there is no conversion from the indexes, it's the initial values
965  // that get used. An alternative would be to only set the pixel
966  // values (no indexes), then set the image as PseudoClass and GM
967  // would create a colormap for us. However, we wouldn't have control
968  // over the order of that colormap. And that's why we set both.
969  Magick::PixelPacket* pix = m_img.getPixels (0, 0, nCols, nRows);
970  Magick::IndexPacket* ind = m_img.getIndexes ();
971  const P* img_fvec = img.fortran_vec ();
972 
973  octave_idx_type GM_idx = 0;
974  for (octave_idx_type column = 0; column < nCols; column++)
975  {
976  for (octave_idx_type row = 0; row < nRows; row++)
977  {
978  ind[GM_idx] = double (*img_fvec);
979  pix[GM_idx] = m_img.colorMap (double (*img_fvec));
980  img_fvec++;
981  GM_idx += nCols;
982  }
983  GM_idx -= nCols * nRows - 1;
984  }
985 
986  // Save changes to underlying image.
987  m_img.syncPixels ();
988  imvec.push_back (m_img);
989  }
990 }
991 
992 static void
993 encode_bool_image (std::vector<Magick::Image>& imvec, const boolNDArray& img)
994 {
995  const octave_idx_type nFrames = img.ndims () < 4 ? 1 : img.dims ()(3);
996  const octave_idx_type nRows = img.rows ();
997  const octave_idx_type nCols = img.columns ();
998 
999  // The initialized image will be black, this is for the other pixels
1000  const Magick::Color white ("white");
1001 
1002  const bool *img_fvec = img.fortran_vec ();
1003  octave_idx_type img_idx = 0;
1004  for (octave_idx_type frame = 0; frame < nFrames; frame++)
1005  {
1006  OCTAVE_QUIT;
1007  // For some reason, we can't set the type to Magick::BilevelType or
1008  // the output image will be black, changing to white has no effect.
1009  // However, this will still work fine and a binary image will be
1010  // saved because we are setting the bitdepth to 1.
1011  Magick::Image m_img = init_enconde_image (nCols, nRows, 1,
1012  Magick::GrayscaleType,
1013  Magick::DirectClass);
1014 
1015  Magick::PixelPacket *pix = m_img.getPixels (0, 0, nCols, nRows);
1016  octave_idx_type GM_idx = 0;
1017  for (octave_idx_type col = 0; col < nCols; col++)
1018  {
1019  for (octave_idx_type row = 0; row < nRows; row++)
1020  {
1021  if (img_fvec[img_idx])
1022  pix[GM_idx] = white;
1023 
1024  img_idx++;
1025  GM_idx += nCols;
1026  }
1027  GM_idx -= nCols * nRows - 1;
1028  }
1029  // Save changes to underlying image.
1030  m_img.syncPixels ();
1031  // While we could not set it to Bilevel at the start, we can do it
1032  // here otherwise some coders won't save it as binary.
1033  m_img.type (Magick::BilevelType);
1034  imvec.push_back (m_img);
1035  }
1036 }
1037 
1038 template <class T>
1039 static void
1040 encode_uint_image (std::vector<Magick::Image>& imvec,
1041  const T& img, const T& alpha)
1042 {
1043  typedef typename T::element_type P;
1044  const octave_idx_type channels = img.ndims () < 3 ? 1 : img.dims ()(2);
1045  const octave_idx_type nFrames = img.ndims () < 4 ? 1 : img.dims ()(3);
1046  const octave_idx_type nRows = img.rows ();
1047  const octave_idx_type nCols = img.columns ();
1048  const octave_idx_type bitdepth = bitdepth_from_class<T> ();
1049 
1050  Magick::ImageType type;
1051  const bool has_alpha = ! alpha.is_empty ();
1052  switch (channels)
1053  {
1054  case 1:
1055  if (has_alpha)
1056  type = Magick::GrayscaleMatteType;
1057  else
1058  type = Magick::GrayscaleType;
1059  break;
1060 
1061  case 3:
1062  if (has_alpha)
1063  type = Magick::TrueColorMatteType;
1064  else
1065  type = Magick::TrueColorType;
1066  break;
1067 
1068  case 4:
1069  if (has_alpha)
1070  type = Magick::ColorSeparationMatteType;
1071  else
1072  type = Magick::ColorSeparationType;
1073  break;
1074 
1075  default:
1076  {
1077  // __imwrite should have already filtered this cases
1078  error ("__magick_write__: wrong size on 3rd dimension");
1079  return;
1080  }
1081  }
1082 
1083  // We will be passing the values as integers with depth as specified
1084  // by QuantumDepth (maximum value specified by MaxRGB). This is independent
1085  // of the actual depth of the image. GM will then convert the values but
1086  // while in memory, it always keeps the values as specified by QuantumDepth.
1087  // From GM documentation:
1088  // Color arguments are must be scaled to fit the Quantum size according to
1089  // the range of MaxRGB
1090  const double divisor = static_cast<double>((uint64_t (1) << bitdepth) - 1)
1091  / MaxRGB;
1092 
1093  const P *img_fvec = img.fortran_vec ();
1094  const P *a_fvec = alpha.fortran_vec ();
1095  switch (type)
1096  {
1097  case Magick::GrayscaleType:
1098  {
1099  for (octave_idx_type frame = 0; frame < nFrames; frame++)
1100  {
1101  OCTAVE_QUIT;
1102  Magick::Image m_img = init_enconde_image (nCols, nRows, bitdepth,
1103  type,
1104  Magick::DirectClass);
1105 
1106  Magick::PixelPacket *pix = m_img.getPixels (0, 0, nCols, nRows);
1107  octave_idx_type GM_idx = 0;
1108  for (octave_idx_type col = 0; col < nCols; col++)
1109  {
1110  for (octave_idx_type row = 0; row < nRows; row++)
1111  {
1112  const double grey = double (*img_fvec) / divisor;
1113  Magick::Color c (grey, grey, grey);
1114  pix[GM_idx] = c;
1115  img_fvec++;
1116  GM_idx += nCols;
1117  }
1118  GM_idx -= nCols * nRows - 1;
1119  }
1120  // Save changes to underlying image.
1121  m_img.syncPixels ();
1122  imvec.push_back (m_img);
1123  }
1124  break;
1125  }
1126 
1127  case Magick::GrayscaleMatteType:
1128  {
1129  for (octave_idx_type frame = 0; frame < nFrames; frame++)
1130  {
1131  OCTAVE_QUIT;
1132  Magick::Image m_img = init_enconde_image (nCols, nRows, bitdepth,
1133  type,
1134  Magick::DirectClass);
1135 
1136  Magick::PixelPacket *pix = m_img.getPixels (0, 0, nCols, nRows);
1137  octave_idx_type GM_idx = 0;
1138  for (octave_idx_type col = 0; col < nCols; col++)
1139  {
1140  for (octave_idx_type row = 0; row < nRows; row++)
1141  {
1142  double grey = double (*img_fvec) / divisor;
1143  Magick::Color c (grey, grey, grey,
1144  MaxRGB - (double (*a_fvec) / divisor));
1145  pix[GM_idx] = c;
1146  img_fvec++;
1147  a_fvec++;
1148  GM_idx += nCols;
1149  }
1150  GM_idx -= nCols * nRows - 1;
1151  }
1152  // Save changes to underlying image.
1153  m_img.syncPixels ();
1154  imvec.push_back (m_img);
1155  }
1156  break;
1157  }
1158 
1159  case Magick::TrueColorType:
1160  {
1161  // The fortran_vec offset for the green and blue channels
1162  const octave_idx_type G_offset = nCols * nRows;
1163  const octave_idx_type B_offset = nCols * nRows * 2;
1164  for (octave_idx_type frame = 0; frame < nFrames; frame++)
1165  {
1166  OCTAVE_QUIT;
1167  Magick::Image m_img = init_enconde_image (nCols, nRows, bitdepth,
1168  type,
1169  Magick::DirectClass);
1170 
1171  Magick::PixelPacket *pix = m_img.getPixels (0, 0, nCols, nRows);
1172  octave_idx_type GM_idx = 0;
1173  for (octave_idx_type col = 0; col < nCols; col++)
1174  {
1175  for (octave_idx_type row = 0; row < nRows; row++)
1176  {
1177  Magick::Color c (double (*img_fvec) / divisor,
1178  double (img_fvec[G_offset]) / divisor,
1179  double (img_fvec[B_offset]) / divisor);
1180  pix[GM_idx] = c;
1181  img_fvec++;
1182  GM_idx += nCols;
1183  }
1184  GM_idx -= nCols * nRows - 1;
1185  }
1186  // Save changes to underlying image.
1187  m_img.syncPixels ();
1188  imvec.push_back (m_img);
1189  img_fvec += B_offset;
1190  }
1191  break;
1192  }
1193 
1194  case Magick::TrueColorMatteType:
1195  {
1196  // The fortran_vec offset for the green and blue channels
1197  const octave_idx_type G_offset = nCols * nRows;
1198  const octave_idx_type B_offset = nCols * nRows * 2;
1199  for (octave_idx_type frame = 0; frame < nFrames; frame++)
1200  {
1201  OCTAVE_QUIT;
1202  Magick::Image m_img = init_enconde_image (nCols, nRows, bitdepth,
1203  type,
1204  Magick::DirectClass);
1205 
1206  Magick::PixelPacket *pix = m_img.getPixels (0, 0, nCols, nRows);
1207  octave_idx_type GM_idx = 0;
1208  for (octave_idx_type col = 0; col < nCols; col++)
1209  {
1210  for (octave_idx_type row = 0; row < nRows; row++)
1211  {
1212  Magick::Color c (double (*img_fvec) / divisor,
1213  double (img_fvec[G_offset]) / divisor,
1214  double (img_fvec[B_offset]) / divisor,
1215  MaxRGB - (double (*a_fvec) / divisor));
1216  pix[GM_idx] = c;
1217  img_fvec++;
1218  a_fvec++;
1219  GM_idx += nCols;
1220  }
1221  GM_idx -= nCols * nRows - 1;
1222  }
1223  // Save changes to underlying image.
1224  m_img.syncPixels ();
1225  imvec.push_back (m_img);
1226  img_fvec += B_offset;
1227  }
1228  break;
1229  }
1230 
1231  case Magick::ColorSeparationType:
1232  {
1233  // The fortran_vec offset for the Magenta, Yellow, and blacK channels
1234  const octave_idx_type M_offset = nCols * nRows;
1235  const octave_idx_type Y_offset = nCols * nRows * 2;
1236  const octave_idx_type K_offset = nCols * nRows * 3;
1237  for (octave_idx_type frame = 0; frame < nFrames; frame++)
1238  {
1239  OCTAVE_QUIT;
1240  Magick::Image m_img = init_enconde_image (nCols, nRows, bitdepth,
1241  type,
1242  Magick::DirectClass);
1243 
1244  Magick::PixelPacket *pix = m_img.getPixels (0, 0, nCols, nRows);
1245  octave_idx_type GM_idx = 0;
1246  for (octave_idx_type col = 0; col < nCols; col++)
1247  {
1248  for (octave_idx_type row = 0; row < nRows; row++)
1249  {
1250  Magick::Color c (double (*img_fvec) / divisor,
1251  double (img_fvec[M_offset]) / divisor,
1252  double (img_fvec[Y_offset]) / divisor,
1253  double (img_fvec[K_offset]) / divisor);
1254  pix[GM_idx] = c;
1255  img_fvec++;
1256  GM_idx += nCols;
1257  }
1258  GM_idx -= nCols * nRows - 1;
1259  }
1260  // Save changes to underlying image.
1261  m_img.syncPixels ();
1262  imvec.push_back (m_img);
1263  img_fvec += K_offset;
1264  }
1265  break;
1266  }
1267 
1268  case Magick::ColorSeparationMatteType:
1269  {
1270  // The fortran_vec offset for the Magenta, Yellow, and blacK channels
1271  const octave_idx_type M_offset = nCols * nRows;
1272  const octave_idx_type Y_offset = nCols * nRows * 2;
1273  const octave_idx_type K_offset = nCols * nRows * 3;
1274  for (octave_idx_type frame = 0; frame < nFrames; frame++)
1275  {
1276  OCTAVE_QUIT;
1277  Magick::Image m_img = init_enconde_image (nCols, nRows, bitdepth,
1278  type,
1279  Magick::DirectClass);
1280 
1281  Magick::PixelPacket *pix = m_img.getPixels (0, 0, nCols, nRows);
1282  Magick::IndexPacket *ind = m_img.getIndexes ();
1283  octave_idx_type GM_idx = 0;
1284  for (octave_idx_type col = 0; col < nCols; col++)
1285  {
1286  for (octave_idx_type row = 0; row < nRows; row++)
1287  {
1288  Magick::Color c (double (*img_fvec) / divisor,
1289  double (img_fvec[M_offset]) / divisor,
1290  double (img_fvec[Y_offset]) / divisor,
1291  double (img_fvec[K_offset]) / divisor);
1292  pix[GM_idx] = c;
1293  ind[GM_idx] = MaxRGB - (double (*a_fvec) / divisor);
1294  img_fvec++;
1295  a_fvec++;
1296  GM_idx += nCols;
1297  }
1298  GM_idx -= nCols * nRows - 1;
1299  }
1300  // Save changes to underlying image.
1301  m_img.syncPixels ();
1302  imvec.push_back (m_img);
1303  img_fvec += K_offset;
1304  }
1305  break;
1306  }
1307 
1308  default:
1309  {
1310  error ("__magick_write__: unrecognized Magick::ImageType");
1311  return;
1312  }
1313  }
1314  return;
1315 }
1316 
1317 // Meant to be shared with both imfinfo and imwrite.
1318 static std::map<octave_idx_type, std::string>
1320 {
1321  // GIF Specifications:
1322  //
1323  // Disposal Method - Indicates the way in which the graphic is to
1324  // be treated after being displayed.
1325  //
1326  // 0 - No disposal specified. The decoder is
1327  // not required to take any action.
1328  // 1 - Do not dispose. The graphic is to be left
1329  // in place.
1330  // 2 - Restore to background color. The area used by the
1331  // graphic must be restored to the background color.
1332  // 3 - Restore to previous. The decoder is required to
1333  // restore the area overwritten by the graphic with
1334  // what was there prior to rendering the graphic.
1335  // 4-7 - To be defined.
1336  static std::map<octave_idx_type, std::string> methods;
1337  if (methods.empty ())
1338  {
1339  methods[0] = "doNotSpecify";
1340  methods[1] = "leaveInPlace";
1341  methods[2] = "restoreBG";
1342  methods[3] = "restorePrevious";
1343  }
1344  return methods;
1345 }
1346 static std::map<std::string, octave_idx_type>
1348 {
1349  static std::map<std::string, octave_idx_type> methods;
1350  if (methods.empty ())
1351  {
1352  methods["donotspecify"] = 0;
1353  methods["leaveinplace"] = 1;
1354  methods["restorebg"] = 2;
1355  methods["restoreprevious"] = 3;
1356  }
1357  return methods;
1358 }
1359 
1360 void static
1361 write_file (const std::string& filename,
1362  const std::string& ext,
1363  std::vector<Magick::Image>& imvec)
1364 {
1365  try
1366  {
1367  Magick::writeImages (imvec.begin (), imvec.end (), ext + ":" + filename);
1368  }
1369  catch (Magick::Warning& w)
1370  {
1371  warning ("Magick++ warning: %s", w.what ());
1372  }
1373  catch (Magick::ErrorCoder& e)
1374  {
1375  warning ("Magick++ coder error: %s", e.what ());
1376  }
1377  catch (Magick::Exception& e)
1378  {
1379  error ("Magick++ exception: %s", e.what ());
1380  error_state = 1;
1381  }
1382 }
1383 
1384 #endif
1385 
1386 DEFUN_DLD (__magick_write__, args, ,
1387  "-*- texinfo -*-\n\
1388 @deftypefn {Loadable Function} {} __magick_write__ (@var{fname}, @var{fmt}, @var{img}, @var{map}, @var{options})\n\
1389 Write image with GraphicsMagick or ImageMagick.\n\
1390 \n\
1391 This is a private internal function not intended for direct use.\n\
1392 Use @code{imwrite} instead.\n\
1393 \n\
1394 @seealso{imfinfo, imformats, imread, imwrite}\n\
1395 @end deftypefn")
1396 {
1397  octave_value_list retval;
1398 
1399 #ifndef HAVE_MAGICK
1400  gripe_disabled_feature ("imwrite", "Image IO");
1401 #else
1402 
1404 
1405  if (args.length () != 5 || ! args(0).is_string () || ! args(1).is_string ())
1406  {
1407  print_usage ();
1408  return retval;
1409  }
1410  const std::string filename = args(0).string_value ();
1411  const std::string ext = args(1).string_value ();
1412 
1413  const octave_scalar_map options = args(4).scalar_map_value ();
1414  if (error_state)
1415  {
1416  error ("__magick_write__: OPTIONS must be a struct");
1417  return retval;
1418  }
1419 
1420  const octave_value img = args(2);
1421  const Matrix cmap = args(3).matrix_value ();
1422  if (error_state)
1423  {
1424  error ("__magick_write__: invalid IMG or MAP");
1425  return retval;
1426  }
1427 
1428  std::vector<Magick::Image> imvec;
1429 
1430  if (cmap.is_empty ())
1431  {
1432  const octave_value alpha = options.getfield ("alpha");
1433  if (img.is_bool_type ())
1434  encode_bool_image (imvec, img.bool_array_value ());
1435  else if (img.is_uint8_type ())
1436  encode_uint_image<uint8NDArray> (imvec, img.uint8_array_value (),
1437  alpha.uint8_array_value ());
1438  else if (img.is_uint16_type ())
1439  encode_uint_image<uint16NDArray> (imvec, img.uint16_array_value (),
1440  alpha.uint16_array_value ());
1441  else if (img.is_uint32_type ())
1442  encode_uint_image<uint32NDArray> (imvec, img.uint32_array_value (),
1443  alpha.uint32_array_value ());
1444  else if (img.is_float_type ())
1445  {
1446  // For image formats that support floating point values, we write
1447  // the actual values. For those who don't, we only use the values
1448  // on the range [0 1] and save integer values.
1449  // But here, even for formats that would support floating point
1450  // values, GM seems unable to do that so we at least make them uint32.
1451  uint32NDArray clip_img;
1452  uint32NDArray clip_alpha;
1453  if (img.is_single_type ())
1454  {
1455  clip_img = img_float2uint<FloatNDArray>
1456  (img.float_array_value ());
1457  clip_alpha = img_float2uint<FloatNDArray>
1458  (alpha.float_array_value ());
1459  }
1460  else
1461  {
1462  clip_img = img_float2uint<NDArray> (img.array_value ());
1463  clip_alpha = img_float2uint<NDArray> (alpha.array_value ());
1464  }
1465  encode_uint_image<uint32NDArray> (imvec, clip_img, clip_alpha);
1466  }
1467  else
1468  {
1469  error ("__magick_write__: image type not supported");
1470  return retval;
1471  }
1472  }
1473  else
1474  {
1475  // We should not get floating point indexed images here because we
1476  // converted them in __imwrite__.m. We should probably do it here
1477  // but it would look much messier.
1478  if (img.is_uint8_type ())
1479  encode_indexed_images<uint8NDArray> (imvec, img.uint8_array_value (),
1480  cmap);
1481  else if (img.is_uint16_type ())
1482  encode_indexed_images<uint16NDArray> (imvec, img.uint16_array_value (),
1483  cmap);
1484  else
1485  {
1486  error ("__magick_write__: indexed image must be uint8, uint16 or float.");
1487  return retval;
1488  }
1489  }
1490  static std::map<std::string, octave_idx_type> disposal_methods
1492 
1493  const octave_idx_type nFrames = imvec.size ();
1494 
1495  const octave_idx_type quality = options.getfield ("quality").int_value ();
1496  const ColumnVector delaytime =
1497  options.getfield ("delaytime").column_vector_value ();
1498  const Array<std::string> disposalmethod =
1499  options.getfield ("disposalmethod").cellstr_value ();
1500  for (octave_idx_type i = 0; i < nFrames; i++)
1501  {
1502  imvec[i].quality (quality);
1503  imvec[i].animationDelay (delaytime(i));
1504  imvec[i].gifDisposeMethod (disposal_methods[disposalmethod(i)]);
1505  }
1506 
1507  // If writemode is set to append, read the image and append to it. Even
1508  // if set to append, make sure that something was read at all.
1509  const std::string writemode = options.getfield ("writemode").string_value ();
1510  if (writemode == "append" && file_stat (filename).exists ())
1511  {
1512  std::vector<Magick::Image> ini_imvec;
1513  read_file (filename, ini_imvec);
1514  if (error_state)
1515  return retval;
1516  if (ini_imvec.size () > 0)
1517  {
1518  ini_imvec.insert (ini_imvec.end (), imvec.begin (), imvec.end ());
1519  ini_imvec.swap (imvec);
1520  }
1521  }
1522 
1523  // FIXME: LoopCount or animationIterations
1524  // How it should work:
1525  //
1526  // This value is only set for the first image in the sequence. Trying
1527  // to set this value with the append mode should have no effect, the
1528  // value used with the first image is the one that counts (that would
1529  // also be Matlab compatible). Thus, the right way to do this would be
1530  // to have an else block on the condition above, and set this only
1531  // when creating a new file. Since Matlab does not interpret a 4D
1532  // matrix as sequence of images to write, its users need to use a for
1533  // loop and set LoopCount only on the first iteration (it actually
1534  // throws warnings otherwise)
1535  //
1536  // Why is this not done the right way:
1537  //
1538  // When GM saves a single image, it discards the value if there is only
1539  // a single image and sets it to "no loop". Since our default is an
1540  // infinite loop, if the user tries to do it the Matlab way (setting
1541  // LoopCount only on the first image) that value will go nowhere.
1542  // See https://sourceforge.net/p/graphicsmagick/bugs/248/
1543  // Because of this, we document to set LoopCount on every iteration
1544  // (in Matlab will cause a lot of warnings), or pass a 4D matrix with
1545  // all frames (won't work in Matlab at all).
1546  // Note that this only needs to be set on the first frame
1547  imvec[0].animationIterations (options.getfield ("loopcount").uint_value ());
1548 
1549  write_file (filename, ext, imvec);
1550  if (error_state)
1551  return retval;
1552 
1553 #endif
1554  return retval;
1555 }
1556 
1557 /*
1558 ## No test needed for internal helper function.
1559 %!assert (1)
1560 */
1561 
1562 // Gets the minimum information from images such as its size and format. Much
1563 // faster than using imfinfo, which slows down a lot since. Note than without
1564 // this, we need to read the image once for imfinfo to set defaults (which is
1565 // done in Octave language), and then again for the actual reading.
1566 DEFUN_DLD (__magick_ping__, args, ,
1567  "-*- texinfo -*-\n\
1568 @deftypefn {Loadable Function} {} __magick_ping__ (@var{fname}, @var{idx})\n\
1569 Ping image information with GraphicsMagick or ImageMagick.\n\
1570 \n\
1571 This is a private internal function not intended for direct use.\n\
1572 \n\
1573 @seealso{imfinfo}\n\
1574 @end deftypefn")
1575 {
1576  octave_value retval;
1577 #ifndef HAVE_MAGICK
1578  gripe_disabled_feature ("imfinfo", "Image IO");
1579 #else
1581 
1582  if (args.length () < 1 || ! args(0).is_string ())
1583  {
1584  print_usage ();
1585  return retval;
1586  }
1587  const std::string filename = args(0).string_value ();
1588  int idx;
1589  if (args.length () > 1)
1590  idx = args(1).int_value () -1;
1591  else
1592  idx = 0;
1593 
1594  Magick::Image img;
1595  img.subImage (idx); // start ping from this image (in case of multi-page)
1596  img.subRange (1); // ping only one of them
1597  try
1598  {
1599  img.ping (filename);
1600  }
1601  catch (Magick::Warning& w)
1602  {
1603  warning ("Magick++ warning: %s", w.what ());
1604  }
1605  catch (Magick::Exception& e)
1606  {
1607  error ("Magick++ exception: %s", e.what ());
1608  return retval;
1609  }
1610 
1611  static const char *fields[] = {"rows", "columns", "format", 0};
1613  ping.setfield ("rows", octave_value (img.rows ()));
1614  ping.setfield ("columns", octave_value (img.columns ()));
1615  ping.setfield ("format", octave_value (img.magick ()));
1616  retval = octave_value (ping);
1617 #endif
1618  return retval;
1619 }
1620 
1621 #ifdef HAVE_MAGICK
1622 static octave_value
1623 magick_to_octave_value (const Magick::CompressionType& magick)
1624 {
1625  switch (magick)
1626  {
1627  case Magick::NoCompression:
1628  return octave_value ("none");
1629  case Magick::BZipCompression:
1630  return octave_value ("bzip");
1631  case Magick::FaxCompression:
1632  return octave_value ("fax3");
1633  case Magick::Group4Compression:
1634  return octave_value ("fax4");
1635  case Magick::JPEGCompression:
1636  return octave_value ("jpeg");
1637  case Magick::LZWCompression:
1638  return octave_value ("lzw");
1639  case Magick::RLECompression:
1640  // This is named "rle" for the HDF, but the same thing is named
1641  // "ccitt" and "PackBits" for binary and non-binary images in TIFF.
1642  return octave_value ("rle");
1643  case Magick::ZipCompression:
1644  return octave_value ("deflate");
1645 
1646  // The following are present only in recent versions of GraphicsMagick.
1647  // At the moment the only use of this would be to have imfinfo report
1648  // the compression method. In the future, someone could implement
1649  // the Compression option for imwrite in which case a macro in
1650  // configure.ac will have to check for their presence of this.
1651  // See bug #39913
1652  // case Magick::LZMACompression:
1653  // return octave_value ("lzma");
1654  // case Magick::JPEG2000Compression:
1655  // return octave_value ("jpeg2000");
1656  // case Magick::JBIG1Compression:
1657  // return octave_value ("jbig1");
1658  // case Magick::JBIG2Compression:
1659  // return octave_value ("jbig2");
1660  default:
1661  return octave_value ("undefined");
1662  }
1663 }
1664 
1665 static octave_value
1666 magick_to_octave_value (const Magick::EndianType& magick)
1667 {
1668  switch (magick)
1669  {
1670  case Magick::LSBEndian:
1671  return octave_value ("little-endian");
1672  case Magick::MSBEndian:
1673  return octave_value ("big-endian");
1674  default:
1675  return octave_value ("undefined");
1676  }
1677 }
1678 
1679 static octave_value
1680 magick_to_octave_value (const Magick::OrientationType& magick)
1681 {
1682  switch (magick)
1683  {
1684  // Values come from the TIFF6 spec
1685  case Magick::TopLeftOrientation:
1686  return octave_value (1);
1687  case Magick::TopRightOrientation:
1688  return octave_value (2);
1689  case Magick::BottomRightOrientation:
1690  return octave_value (3);
1691  case Magick::BottomLeftOrientation:
1692  return octave_value (4);
1693  case Magick::LeftTopOrientation:
1694  return octave_value (5);
1695  case Magick::RightTopOrientation:
1696  return octave_value (6);
1697  case Magick::RightBottomOrientation:
1698  return octave_value (7);
1699  case Magick::LeftBottomOrientation:
1700  return octave_value (8);
1701  default:
1702  return octave_value (1);
1703  }
1704 }
1705 
1706 static octave_value
1707 magick_to_octave_value (const Magick::ResolutionType& magick)
1708 {
1709  switch (magick)
1710  {
1711  case Magick::PixelsPerInchResolution:
1712  return octave_value ("Inch");
1713  case Magick::PixelsPerCentimeterResolution:
1714  return octave_value ("Centimeter");
1715  default:
1716  return octave_value ("undefined");
1717  }
1718 }
1719 
1720 static bool
1721 is_valid_exif (const std::string& val)
1722 {
1723  // Sometimes GM will return the string "unknown" instead of empty
1724  // for an empty value.
1725  return (! val.empty () && val != "unknown");
1726 }
1727 
1728 static void
1729 fill_exif (octave_scalar_map& map, Magick::Image& img,
1730  const std::string& key)
1731 {
1732  const std::string attr = img.attribute ("EXIF:" + key);
1733  if (is_valid_exif (attr))
1734  map.setfield (key, octave_value (attr));
1735  return;
1736 }
1737 
1738 static void
1739 fill_exif_ints (octave_scalar_map& map, Magick::Image& img,
1740  const std::string& key)
1741 {
1742  const std::string attr = img.attribute ("EXIF:" + key);
1743  if (is_valid_exif (attr))
1744  {
1745  // string of the type "float,float,float....."
1746  float number;
1747  ColumnVector values (std::count (attr.begin (), attr.end (), ',') +1);
1748  std::string sub;
1749  std::istringstream sstream (attr);
1750  octave_idx_type n = 0;
1751  while (std::getline (sstream, sub, char (',')))
1752  {
1753  sscanf (sub.c_str (), "%f", &number);
1754  values(n++) = number;
1755  }
1756  map.setfield (key, octave_value (values));
1757  }
1758  return;
1759 }
1760 
1761 static void
1762 fill_exif_floats (octave_scalar_map& map, Magick::Image& img,
1763  const std::string& key)
1764 {
1765  const std::string attr = img.attribute ("EXIF:" + key);
1766  if (is_valid_exif (attr))
1767  {
1768  // string of the type "int/int,int/int,int/int....."
1769  int numerator;
1770  int denominator;
1771  ColumnVector values (std::count (attr.begin (), attr.end (), ',') +1);
1772  std::string sub;
1773  std::istringstream sstream (attr);
1774  octave_idx_type n = 0;
1775  while (std::getline (sstream, sub, ','))
1776  {
1777  sscanf (sub.c_str (), "%i/%i", &numerator, &denominator);
1778  values(n++) = double (numerator) / double (denominator);
1779  }
1780  map.setfield (key, octave_value (values));
1781  }
1782  return;
1783 }
1784 
1785 #endif
1786 
1787 DEFUN_DLD (__magick_finfo__, args, ,
1788  "-*- texinfo -*-\n\
1789 @deftypefn {Loadable Function} {} __magick_finfo__ (@var{fname})\n\
1790 Read image information with GraphicsMagick or ImageMagick.\n\
1791 \n\
1792 This is a private internal function not intended for direct use.\n\
1793 Use @code{imfinfo} instead.\n\
1794 \n\
1795 @seealso{imfinfo, imformats, imread, imwrite}\n\
1796 @end deftypefn")
1797 {
1798  octave_value retval;
1799 
1800 #ifndef HAVE_MAGICK
1801  gripe_disabled_feature ("imfinfo", "Image IO");
1802 #else
1804 
1805  if (args.length () < 1 || ! args(0).is_string ())
1806  {
1807  print_usage ();
1808  return retval;
1809  }
1810  const std::string filename = args(0).string_value ();
1811 
1812  std::vector<Magick::Image> imvec;
1813  read_file (filename, imvec);
1814  if (error_state)
1815  return retval;
1816  const octave_idx_type nFrames = imvec.size ();
1817  const std::string format = imvec[0].magick ();
1818 
1819  // Here's how this function works. We need to return a struct array, one
1820  // struct for each image in the file (remember, there are image
1821  // that allow for multiple images in the same file). Now, Matlab seems
1822  // to have format specific code so the fields on the struct are different
1823  // for each format. It only has a small subset that is common to all
1824  // of them, the others are undocumented. Because we try to abstract from
1825  // the formats we always return the same list of fields (note that with
1826  // GM we support more than 88 formats. That's way more than Matlab, and
1827  // I don't want to write specific code for each of them).
1828  //
1829  // So what we do is we create an octave_scalar_map, fill it with the
1830  // information for that image, and then insert it into an octave_map.
1831  // Because in the same file, different images may have values for
1832  // different fields, we can't create a field only if there's a value.
1833  // Bad things happen if we merge octave_scalar_maps with different
1834  // fields from the others (suppose for example a TIFF file with 4 images,
1835  // where only the third image has a colormap.
1836 
1837  static const char *fields[] =
1838  {
1839  // These are fields that must always appear for Matlab.
1840  "Filename",
1841  "FileModDate",
1842  "FileSize",
1843  "Format",
1844  "FormatVersion",
1845  "Width",
1846  "Height",
1847  "BitDepth",
1848  "ColorType",
1849 
1850  // These are format specific or not existent in Matlab. The most
1851  // annoying thing is that Matlab may have different names for the
1852  // same thing in different formats.
1853  "DelayTime",
1854  "DisposalMethod",
1855  "LoopCount",
1856  "ByteOrder",
1857  "Gamma",
1858  "Chromaticities",
1859  "Comment",
1860  "Quality",
1861  "Compression", // same as CompressionType
1862  "Colormap", // same as ColorTable (in PNG)
1863  "Orientation",
1864  "ResolutionUnit",
1865  "XResolution",
1866  "YResolution",
1867  "Software", // sometimes is an Exif tag
1868  "Make", // actually an Exif tag
1869  "Model", // actually an Exif tag
1870  "DateTime", // actually an Exif tag
1871  "ImageDescription", // actually an Exif tag
1872  "Artist", // actually an Exif tag
1873  "Copyright", // actually an Exif tag
1874  "DigitalCamera",
1875  "GPSInfo",
1876  // Notes for the future: GM allows one to get many attributes, and even has
1877  // attribute() to obtain arbitrary ones, that may exist in only some
1878  // cases. The following is a list of some methods and into what possible
1879  // Matlab compatible values they may be converted.
1880  //
1881  // colorSpace() -> PhotometricInterpretation
1882  // backgroundColor() -> BackgroundColor
1883  // interlaceType() -> Interlaced, InterlaceType, and PlanarConfiguration
1884  // label() -> Title
1885  0
1886  };
1887 
1888  // The one we will return at the end
1889  octave_map info (dim_vector (nFrames, 1), string_vector (fields));
1890 
1891  // Some of the fields in the struct are about file information and will be
1892  // the same for all images in the file. So we create a template, fill in
1893  // those values, and make a copy of the template for each image.
1894  octave_scalar_map template_info = (string_vector (fields));
1895 
1896  template_info.setfield ("Format", octave_value (format));
1897  // We can't actually get FormatVersion but even Matlab sometimes can't.
1898  template_info.setfield ("FormatVersion", octave_value (""));
1899 
1900  const file_stat fs (filename);
1901  if (fs)
1902  {
1903  const octave_localtime mtime (fs.mtime ());
1904  const std::string filetime = mtime.strftime ("%e-%b-%Y %H:%M:%S");
1905  template_info.setfield ("Filename", octave_value (filename));
1906  template_info.setfield ("FileModDate", octave_value (filetime));
1907  template_info.setfield ("FileSize", octave_value (fs.size ()));
1908  }
1909  else
1910  {
1911  error ("imfinfo: error reading '%s': %s",
1912  filename.c_str (), fs.error ().c_str ());
1913  return retval;
1914  }
1915 
1916  for (octave_idx_type frame = 0; frame < nFrames; frame++)
1917  {
1918  OCTAVE_QUIT;
1919  octave_scalar_map info_frame (template_info);
1920  const Magick::Image img = imvec[frame];
1921 
1922  info_frame.setfield ("Width", octave_value (img.columns ()));
1923  info_frame.setfield ("Height", octave_value (img.rows ()));
1924  info_frame.setfield ("BitDepth",
1925  octave_value (get_depth (const_cast<Magick::Image&> (img))));
1926 
1927  // Stuff related to colormap, image class and type
1928  // Because GM is too smart for us... Read the comments in is_indexed()
1929  {
1930  std::string color_type;
1931  Matrix cmap;
1932  if (is_indexed (img))
1933  {
1934  color_type = "indexed";
1935  cmap =
1936  read_maps (const_cast<Magick::Image&> (img))(0).matrix_value ();
1937  }
1938  else
1939  {
1940  switch (img.type ())
1941  {
1942  case Magick::BilevelType:
1943  case Magick::GrayscaleType:
1944  case Magick::GrayscaleMatteType:
1945  color_type = "grayscale";
1946  break;
1947 
1948  case Magick::TrueColorType:
1949  case Magick::TrueColorMatteType:
1950  color_type = "truecolor";
1951  break;
1952 
1953  case Magick::PaletteType:
1954  case Magick::PaletteMatteType:
1955  // we should never get here or is_indexed needs to be fixed
1956  color_type = "indexed";
1957  break;
1958 
1959  case Magick::ColorSeparationType:
1960  case Magick::ColorSeparationMatteType:
1961  color_type = "CMYK";
1962  break;
1963 
1964  default:
1965  color_type = "undefined";
1966  }
1967  }
1968  info_frame.setfield ("ColorType", octave_value (color_type));
1969  info_frame.setfield ("Colormap", octave_value (cmap));
1970  }
1971 
1972  {
1973  // Not all images have chroma values. In such cases, they'll
1974  // be all zeros. So rather than send a matrix of zeros, we will
1975  // check for that, and send an empty vector instead.
1976  RowVector chromaticities (8);
1977  double* chroma_fvec = chromaticities.fortran_vec ();
1978  img.chromaWhitePoint (&chroma_fvec[0], &chroma_fvec[1]);
1979  img.chromaRedPrimary (&chroma_fvec[2], &chroma_fvec[3]);
1980  img.chromaGreenPrimary (&chroma_fvec[4], &chroma_fvec[5]);
1981  img.chromaBluePrimary (&chroma_fvec[6], &chroma_fvec[7]);
1982  if (chromaticities.nnz () == 0)
1983  chromaticities = RowVector (0);
1984  info_frame.setfield ("Chromaticities", octave_value (chromaticities));
1985  }
1986 
1987  info_frame.setfield ("Gamma", octave_value (img.gamma ()));
1988  info_frame.setfield ("XResolution", octave_value (img.xResolution ()));
1989  info_frame.setfield ("YResolution", octave_value (img.yResolution ()));
1990  info_frame.setfield ("DelayTime", octave_value (img.animationDelay ()));
1991  info_frame.setfield ("LoopCount",
1992  octave_value (img.animationIterations ()));
1993  info_frame.setfield ("Quality", octave_value (img.quality ()));
1994  info_frame.setfield ("Comment", octave_value (img.comment ()));
1995 
1996  info_frame.setfield ("Compression",
1997  magick_to_octave_value (img.compressType ()));
1998  info_frame.setfield ("Orientation",
1999  magick_to_octave_value (img.orientation ()));
2000  info_frame.setfield ("ResolutionUnit",
2001  magick_to_octave_value (img.resolutionUnits ()));
2002  info_frame.setfield ("ByteOrder",
2003  magick_to_octave_value (img.endian ()));
2004 
2005  // It is not possible to know if there's an Exif field so we just
2006  // check for the Exif Version value. If it does exists, then we
2007  // bother about looking for specific fields.
2008  {
2009  Magick::Image& cimg = const_cast<Magick::Image&> (img);
2010 
2011  // These will be in Exif tags but must appear as fields in the
2012  // base struct array, not as another struct in one of its fields.
2013  // This is likely because they belong to the Baseline TIFF specs
2014  // and may appear out of the Exif tag. So first we check if it
2015  // exists outside the Exif tag.
2016  // See Section 4.6.4, table 4, page 28 of Exif specs version 2.3
2017  // (CIPA DC- 008-Translation- 2010)
2018  static const char *base_exif_str_fields[] =
2019  {
2020  "DateTime",
2021  "ImageDescription",
2022  "Make",
2023  "Model",
2024  "Software",
2025  "Artist",
2026  "Copyright",
2027  0,
2028  };
2029  static const string_vector base_exif_str (base_exif_str_fields);
2030  static const octave_idx_type n_base_exif_str = base_exif_str.numel ();
2031  for (octave_idx_type field = 0; field < n_base_exif_str; field++)
2032  {
2033  info_frame.setfield (base_exif_str[field],
2034  octave_value (cimg.attribute (base_exif_str[field])));
2035  fill_exif (info_frame, cimg, base_exif_str[field]);
2036  }
2037 
2038  octave_scalar_map camera;
2039  octave_scalar_map gps;
2040  if (! cimg.attribute ("EXIF:ExifVersion").empty ())
2041  {
2042  // See Section 4.6.5, table 7 and 8, over pages page 42 to 43
2043  // of Exif specs version 2.3 (CIPA DC- 008-Translation- 2010)
2044 
2045  // Listed on the Exif specs as being of type ASCII.
2046  static const char *exif_str_fields[] =
2047  {
2048  "RelatedSoundFile",
2049  "DateTimeOriginal",
2050  "DateTimeDigitized",
2051  "SubSecTime",
2052  "DateTimeOriginal",
2053  "SubSecTimeOriginal",
2054  "SubSecTimeDigitized",
2055  "ImageUniqueID",
2056  "CameraOwnerName",
2057  "BodySerialNumber",
2058  "LensMake",
2059  "LensModel",
2060  "LensSerialNumber",
2061  "SpectralSensitivity",
2062  // These last two are of type undefined but most likely will
2063  // be strings. Even if they're not GM returns a string anyway.
2064  "UserComment",
2065  "MakerComment",
2066  0
2067  };
2068  static const string_vector exif_str (exif_str_fields);
2069  static const octave_idx_type n_exif_str = exif_str.numel ();
2070  for (octave_idx_type field = 0; field < n_exif_str; field++)
2071  fill_exif (camera, cimg, exif_str[field]);
2072 
2073  // Listed on the Exif specs as being of type SHORT or LONG.
2074  static const char *exif_int_fields[] =
2075  {
2076  "ColorSpace",
2077  "ExifImageWidth", // PixelXDimension (CPixelXDimension in Matlab)
2078  "ExifImageHeight", // PixelYDimension (CPixelYDimension in Matlab)
2079  "PhotographicSensitivity",
2080  "StandardOutputSensitivity",
2081  "RecommendedExposureIndex",
2082  "ISOSpeed",
2083  "ISOSpeedLatitudeyyy",
2084  "ISOSpeedLatitudezzz",
2085  "FocalPlaneResolutionUnit",
2086  "FocalLengthIn35mmFilm",
2087  // Listed as SHORT or LONG but with more than 1 count.
2088  "SubjectArea",
2089  "SubjectLocation",
2090  // While the following are an integer, their value have a meaning
2091  // that must be represented as a string for Matlab compatibility.
2092  // For example, a 3 on ExposureProgram, would return
2093  // "Aperture priority" as defined on the Exif specs.
2094  "ExposureProgram",
2095  "SensitivityType",
2096  "MeteringMode",
2097  "LightSource",
2098  "Flash",
2099  "SensingMethod",
2100  "FileSource",
2101  "CustomRendered",
2102  "ExposureMode",
2103  "WhiteBalance",
2104  "SceneCaptureType",
2105  "GainControl",
2106  "Contrast",
2107  "Saturation",
2108  "Sharpness",
2109  "SubjectDistanceRange",
2110  0
2111  };
2112  static const string_vector exif_int (exif_int_fields);
2113  static const octave_idx_type n_exif_int = exif_int.numel ();
2114  for (octave_idx_type field = 0; field < n_exif_int; field++)
2115  fill_exif_ints (camera, cimg, exif_int[field]);
2116 
2117  // Listed as RATIONAL or SRATIONAL
2118  static const char *exif_float_fields[] =
2119  {
2120  "Gamma",
2121  "CompressedBitsPerPixel",
2122  "ExposureTime",
2123  "FNumber",
2124  "ShutterSpeedValue", // SRATIONAL
2125  "ApertureValue",
2126  "BrightnessValue", // SRATIONAL
2127  "ExposureBiasValue", // SRATIONAL
2128  "MaxApertureValue",
2129  "SubjectDistance",
2130  "FocalLength",
2131  "FlashEnergy",
2132  "FocalPlaneXResolution",
2133  "FocalPlaneYResolution",
2134  "ExposureIndex",
2135  "DigitalZoomRatio",
2136  // Listed as RATIONAL or SRATIONAL with more than 1 count.
2137  "LensSpecification",
2138  0
2139  };
2140  static const string_vector exif_float (exif_float_fields);
2141  static const octave_idx_type n_exif_float = exif_float.numel ();
2142  for (octave_idx_type field = 0; field < n_exif_float; field++)
2143  fill_exif_floats (camera, cimg, exif_float[field]);
2144 
2145  // Inside a Exif field, it is possible that there is also a
2146  // GPS field. This is not the same as ExifVersion but seems
2147  // to be how we have to check for it.
2148  if (cimg.attribute ("EXIF:GPSInfo") != "unknown")
2149  {
2150  // The story here is the same as with Exif.
2151  // See Section 4.6.6, table 15 on page 68 of Exif specs
2152  // version 2.3 (CIPA DC- 008-Translation- 2010)
2153 
2154  static const char *gps_str_fields[] =
2155  {
2156  "GPSLatitudeRef",
2157  "GPSLongitudeRef",
2158  "GPSAltitudeRef",
2159  "GPSSatellites",
2160  "GPSStatus",
2161  "GPSMeasureMode",
2162  "GPSSpeedRef",
2163  "GPSTrackRef",
2164  "GPSImgDirectionRef",
2165  "GPSMapDatum",
2166  "GPSDestLatitudeRef",
2167  "GPSDestLongitudeRef",
2168  "GPSDestBearingRef",
2169  "GPSDestDistanceRef",
2170  "GPSDateStamp",
2171  0
2172  };
2173  static const string_vector gps_str (gps_str_fields);
2174  static const octave_idx_type n_gps_str = gps_str.numel ();
2175  for (octave_idx_type field = 0; field < n_gps_str; field++)
2176  fill_exif (gps, cimg, gps_str[field]);
2177 
2178  static const char *gps_int_fields[] =
2179  {
2180  "GPSDifferential",
2181  0
2182  };
2183  static const string_vector gps_int (gps_int_fields);
2184  static const octave_idx_type n_gps_int = gps_int.numel ();
2185  for (octave_idx_type field = 0; field < n_gps_int; field++)
2186  fill_exif_ints (gps, cimg, gps_int[field]);
2187 
2188  static const char *gps_float_fields[] =
2189  {
2190  "GPSAltitude",
2191  "GPSDOP",
2192  "GPSSpeed",
2193  "GPSTrack",
2194  "GPSImgDirection",
2195  "GPSDestBearing",
2196  "GPSDestDistance",
2197  "GPSHPositioningError",
2198  // Listed as RATIONAL or SRATIONAL with more than 1 count.
2199  "GPSLatitude",
2200  "GPSLongitude",
2201  "GPSTimeStamp",
2202  "GPSDestLatitude",
2203  "GPSDestLongitude",
2204  0
2205  };
2206  static const string_vector gps_float (gps_float_fields);
2207  static const octave_idx_type n_gps_float = gps_float.numel ();
2208  for (octave_idx_type field = 0; field < n_gps_float; field++)
2209  fill_exif_floats (gps, cimg, gps_float[field]);
2210 
2211  }
2212  }
2213  info_frame.setfield ("DigitalCamera", octave_value (camera));
2214  info_frame.setfield ("GPSInfo", octave_value (gps));
2215  }
2216 
2217  info.fast_elem_insert (frame, info_frame);
2218  }
2219 
2220  if (format == "GIF")
2221  {
2222  static std::map<octave_idx_type, std::string> disposal_methods
2223  = init_disposal_methods ();
2224  string_vector methods (nFrames);
2225  for (octave_idx_type frame = 0; frame < nFrames; frame++)
2226  methods[frame] = disposal_methods[imvec[frame].gifDisposeMethod ()];
2227  info.setfield ("DisposalMethod", Cell (methods));
2228  }
2229  else
2230  info.setfield ("DisposalMethod",
2231  Cell (dim_vector (nFrames, 1), octave_value ("")));
2232 
2233  retval = octave_value (info);
2234 #endif
2235  return retval;
2236 }
2237 
2238 /*
2239 ## No test needed for internal helper function.
2240 %!assert (1)
2241 */
2242 
2243 DEFUN_DLD (__magick_formats__, args, ,
2244  "-*- texinfo -*-\n\
2245 @deftypefn {Loadable Function} {} __magick_imformats__ (@var{formats})\n\
2246 Fill formats info with GraphicsMagick CoderInfo.\n\
2247 \n\
2248 @seealso{imfinfo, imformats, imread, imwrite}\n\
2249 @end deftypefn")
2250 {
2251  octave_value retval;
2252 #ifndef HAVE_MAGICK
2253  gripe_disabled_feature ("imformats", "Image IO");
2254 #else
2255  if (args.length () != 1 || ! args(0).is_map ())
2256  {
2257  print_usage ();
2258  return retval;
2259  }
2260  octave_map formats = args(0).map_value ();
2261 
2263  for (octave_idx_type idx = 0; idx < formats.numel (); idx++)
2264  {
2265  try
2266  {
2267  octave_scalar_map fmt = formats.checkelem (idx);
2268  Magick::CoderInfo coder (fmt.getfield ("coder").string_value ());
2269 
2270  fmt.setfield ("description", octave_value (coder.description ()));
2271  fmt.setfield ("multipage", coder.isMultiFrame () ? true : false);
2272  // default for read and write is a function handle. If we can't
2273  // read or write them, them set it to an empty value
2274  if (! coder.isReadable ())
2275  fmt.setfield ("read", Matrix ());
2276  if (! coder.isWritable ())
2277  fmt.setfield ("write", Matrix ());
2278  formats.fast_elem_insert (idx, fmt);
2279  }
2280  catch (Magick::Exception& e)
2281  {
2282  // Exception here are missing formats. So we remove the format
2283  // from the structure and reduce idx.
2284  formats.delete_elements (idx);
2285  idx--;
2286  }
2287  }
2288  retval = formats;
2289 #endif
2290  return retval;
2291 }
2292 
2293 /*
2294 ## No test needed for internal helper function.
2295 %!assert (1)
2296 */
uint8NDArray uint8_array_value(void) const
Definition: ov.h:882
void warning_with_id(const char *id, const char *fmt,...)
Definition: error.cc:696
ColumnVector column_vector_value(bool frc_str_conv=false, bool frc_vec_conv=false) const
Definition: ov.cc:1658
bool is_empty(void) const
Definition: Array.h:472
bool is_range(void) const
Definition: ov.h:571
std::string error(void) const
Definition: file-stat.h:136
Definition: Cell.h:35
void delete_elements(const idx_vector &i)
Definition: oct-map.cc:1179
unsigned int uint_value(bool req_int=false, bool frc_str_conv=false) const
Definition: ov.h:734
bool is_uint16_type(void) const
Definition: ov.h:634
int ndims(void) const
Definition: Array.h:487
OCTINTERP_API void print_usage(void)
Definition: defun.cc:51
octave_idx_type numel(void) const
Number of elements in the array.
Definition: Array.h:275
static Range get_region_range(const octave_value &region)
std::string strftime(const std::string &fmt) const
Definition: oct-time.cc:144
octave_idx_type nelem(void) const
Definition: Range.h:64
static bool is_indexed(const Magick::Image &img)
bool is_scalar_type(void) const
Definition: ov.h:657
int int_value(bool req_int=false, bool frc_str_conv=false) const
Definition: ov.h:730
Definition: Range.h:31
static std::map< std::string, octave_idx_type > init_reverse_disposal_methods()
void error(const char *fmt,...)
Definition: error.cc:476
static void read_file(const std::string &filename, std::vector< Magick::Image > &imvec)
void setfield(const std::string &key, const octave_value &val)
Definition: oct-map.cc:171
static octave_idx_type bitdepth_from_class()
static void write_file(const std::string &filename, const std::string &ext, std::vector< Magick::Image > &imvec)
static void maybe_initialize_magick(void)
octave_idx_type numel(void) const
Definition: oct-map.h:372
double scalar_value(bool frc_str_conv=false) const
Definition: ov.h:765
static void encode_bool_image(std::vector< Magick::Image > &imvec, const boolNDArray &img)
boolNDArray bool_array_value(bool warn=false) const
Definition: ov.h:811
octave_idx_type rows(void) const
Definition: Array.h:313
static octave_value magick_to_octave_value(const Magick::CompressionType &magick)
octave_idx_type nelem(void) const
Number of elements in the array.
Definition: Array.h:271
static octave_value_list read_indexed_images(const std::vector< Magick::Image > &imvec, const Array< octave_idx_type > &frameidx, const octave_idx_type &nargout, const octave_scalar_map &options)
Cell cell_value(void) const
Definition: ov.cc:1566
bool is_float_type(void) const
Definition: ov.h:614
Array< std::string > cellstr_value(void) const
Definition: ov.h:900
const dim_vector & dims(void) const
Return a const-reference so that dims ()(i) works efficiently.
Definition: Array.h:337
static std::string get_program_invocation_name(void)
Definition: oct-env.cc:167
static bool is_valid_exif(const std::string &val)
FloatNDArray float_array_value(bool frc_str_conv=false) const
Definition: ov.h:782
double max(void) const
Definition: Range.cc:215
bool is_bool_type(void) const
Definition: ov.h:645
std::string string_value(bool force=false) const
Definition: ov.h:897
Range range_value(void) const
Definition: ov.h:903
std::complex< double > w(std::complex< double > z, double relerr=0)
void gripe_disabled_feature(const std::string &func, const std::string &feature, const std::string &pkg)
Definition: gripes.cc:242
bool is_string(void) const
Definition: ov.h:562
double inc(void) const
Definition: Range.h:63
int error_state
Definition: error.cc:101
bool fast_elem_insert(octave_idx_type n, const octave_scalar_map &rhs)
Definition: oct-map.cc:400
void resize(const dim_vector &dv, const T &rfv)
Definition: Array.cc:1033
static octave_int< T > max(void)
Definition: oct-inttypes.h:932
static void encode_indexed_images(std::vector< Magick::Image > &imvec, const T &img, const Matrix &cmap)
Definition: dMatrix.h:35
static void fill_exif_ints(octave_scalar_map &map, Magick::Image &img, const std::string &key)
static void fill_exif(octave_scalar_map &map, Magick::Image &img, const std::string &key)
octave_idx_type nnz(void) const
Count nonzero elements.
Definition: Array.cc:2231
static octave_value_list read_maps(Magick::Image &img)
octave_value_list read_images(std::vector< Magick::Image > &imvec, const Array< octave_idx_type > &frameidx, const octave_idx_type &nargout, const octave_scalar_map &options)
static Magick::Image init_enconde_image(const octave_idx_type &nCols, const octave_idx_type &nRows, const octave_idx_type &bitdepth, const Magick::ImageType &type, const Magick::ClassType &klass)
void setfield(const std::string &key, const Cell &val)
Definition: oct-map.cc:265
void warning(const char *fmt,...)
Definition: error.cc:681
charNDArray max(char d, const charNDArray &m)
Definition: chNDArray.cc:233
octave_idx_type length(void) const
Number of elements in the array.
Definition: Array.h:267
octave_time mtime(void) const
Definition: file-stat.h:118
NDArray array_value(bool frc_str_conv=false) const
Definition: ov.h:779
static std::map< std::string, octave_idx_type > calculate_region(const octave_scalar_map &options)
static octave_idx_type get_depth(Magick::Image &img)
bool is_uint8_type(void) const
Definition: ov.h:631
static uint32NDArray img_float2uint(const T &img)
static std::map< octave_idx_type, std::string > init_disposal_methods()
static void encode_uint_image(std::vector< Magick::Image > &imvec, const T &img, const T &alpha)
double base(void) const
Definition: Range.h:61
octave_scalar_map checkelem(octave_idx_type n) const
Definition: oct-map.cc:357
octave_value getfield(const std::string &key) const
Definition: oct-map.cc:164
Array< int > int_vector_value(bool req_int=false, bool frc_str_conv=false, bool frc_vec_conv=false) const
Definition: ov.cc:1717
#define DEFUN_DLD(name, args_name, nargout_name, doc)
Definition: defun-dld.h:59
const T * fortran_vec(void) const
Definition: Array.h:481
bool is_single_type(void) const
Definition: ov.h:611
static void fill_exif_floats(octave_scalar_map &map, Magick::Image &img, const std::string &key)
#define OCTAVE_QUIT
Definition: quit.h:130
off_t size(void) const
Definition: file-stat.h:115
bool is_uint32_type(void) const
Definition: ov.h:637
uint32NDArray uint32_array_value(void) const
Definition: ov.h:888
uint16NDArray uint16_array_value(void) const
Definition: ov.h:885
octave_idx_type columns(void) const
Definition: Array.h:322
return octave_value(v1.char_array_value().concat(v2.char_array_value(), ra_idx),((a1.is_sq_string()||a2.is_sq_string())? '\'': '"'))