FFmpeg
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
j2k.c
Go to the documentation of this file.
1 /*
2  * JPEG 2000 encoder and decoder common functions
3  * Copyright (c) 2007 Kamil Nowosad
4  * Copyright (c) 2013 Nicolas Bertrand <nicoinattendu@gmail.com>
5  *
6  * This file is part of FFmpeg.
7  *
8  * FFmpeg is free software; you can redistribute it and/or
9  * modify it under the terms of the GNU Lesser General Public
10  * License as published by the Free Software Foundation; either
11  * version 2.1 of the License, or (at your option) any later version.
12  *
13  * FFmpeg is distributed in the hope that it will be useful,
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16  * Lesser General Public License for more details.
17  *
18  * You should have received a copy of the GNU Lesser General Public
19  * License along with FFmpeg; if not, write to the Free Software
20  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
21  */
22 
23 /**
24  * JPEG 2000 image encoder and decoder common functions
25  * @file
26  * @author Kamil Nowosad
27  */
28 
29 #include "libavutil/avassert.h"
30 #include "libavutil/common.h"
31 #include "libavutil/mem.h"
32 #include "avcodec.h"
33 #include "j2k.h"
34 
35 #define SHL(a, n) ((n) >= 0 ? (a) << (n) : (a) >> -(n))
36 
37 /* tag tree routines */
38 
39 /* allocate the memory for tag tree */
40 static int tag_tree_size(int w, int h)
41 {
42  int res = 0;
43  while (w > 1 || h > 1) {
44  res += w * h;
45  w = (w + 1) >> 1;
46  h = (h + 1) >> 1;
47  }
48  return res + 1;
49 }
50 
52 {
53  int pw = w, ph = h;
54  Jpeg2000TgtNode *res, *t, *t2;
55  int32_t tt_size;
56 
57  tt_size = tag_tree_size(w, h);
58 
59  t = res = av_mallocz_array(tt_size, sizeof(*t));
60  if (!res)
61  return NULL;
62 
63  while (w > 1 || h > 1) {
64  int i, j;
65  pw = w;
66  ph = h;
67 
68  w = (w + 1) >> 1;
69  h = (h + 1) >> 1;
70  t2 = t + pw * ph;
71 
72  for (i = 0; i < ph; i++)
73  for (j = 0; j < pw; j++)
74  t[i * pw + j].parent = &t2[(i >> 1) * w + (j >> 1)];
75 
76  t = t2;
77  }
78  t[0].parent = NULL;
79  return res;
80 }
81 
82 static void tag_tree_zero(Jpeg2000TgtNode *t, int w, int h)
83 {
84  int i, siz = tag_tree_size(w, h);
85 
86  for (i = 0; i < siz; i++) {
87  t[i].val = 0;
88  t[i].vis = 0;
89  }
90 }
91 
92 static int getsigctxno(int flag, int bandno)
93 {
94  int h, v, d;
95 
96  h = ((flag & JPEG2000_T1_SIG_E) ? 1 : 0) +
97  ((flag & JPEG2000_T1_SIG_W) ? 1 : 0);
98  v = ((flag & JPEG2000_T1_SIG_N) ? 1 : 0) +
99  ((flag & JPEG2000_T1_SIG_S) ? 1 : 0);
100  d = ((flag & JPEG2000_T1_SIG_NE) ? 1 : 0) +
101  ((flag & JPEG2000_T1_SIG_NW) ? 1 : 0) +
102  ((flag & JPEG2000_T1_SIG_SE) ? 1 : 0) +
103  ((flag & JPEG2000_T1_SIG_SW) ? 1 : 0);
104 
105  if (bandno < 3) {
106  if (bandno == 1)
107  FFSWAP(int, h, v);
108  if (h == 2) return 8;
109  if (h == 1) {
110  if (v >= 1) return 7;
111  if (d >= 1) return 6;
112  return 5;
113  }
114  if (v == 2) return 4;
115  if (v == 1) return 3;
116  if (d >= 2) return 2;
117  if (d == 1) return 1;
118  } else {
119  if (d >= 3) return 8;
120  if (d == 2) {
121  if (h+v >= 1) return 7;
122  return 6;
123  }
124  if (d == 1) {
125  if (h+v >= 2) return 5;
126  if (h+v == 1) return 4;
127  return 3;
128  }
129  if (h+v >= 2) return 2;
130  if (h+v == 1) return 1;
131  }
132  return 0;
133 }
134 
135 
136 static const int contribtab[3][3] = { { 0, -1, 1 }, { -1, -1, 0 }, { 1, 0, 1 } };
137 static const int ctxlbltab[3][3] = { { 13, 12, 11 }, { 10, 9, 10 }, { 11, 12, 13 } };
138 static const int xorbittab[3][3] = { { 1, 1, 1 }, { 1, 0, 0 }, { 0, 0, 0 } };
139 
140 static int getsgnctxno(int flag, uint8_t *xorbit)
141 {
142  int vcontrib, hcontrib;
143 
144  hcontrib = contribtab[flag & JPEG2000_T1_SIG_E ? flag & JPEG2000_T1_SGN_E ? 1 : 2 : 0]
145  [flag & JPEG2000_T1_SIG_W ? flag & JPEG2000_T1_SGN_W ? 1 : 2 : 0] + 1;
146  vcontrib = contribtab[flag & JPEG2000_T1_SIG_S ? flag & JPEG2000_T1_SGN_S ? 1 : 2 : 0]
147  [flag & JPEG2000_T1_SIG_N ? flag & JPEG2000_T1_SGN_N ? 1 : 2 : 0] + 1;
148  *xorbit = xorbittab[hcontrib][vcontrib];
149 
150  return ctxlbltab[hcontrib][vcontrib];
151 }
152 
154  int negative)
155 {
156  x++;
157  y++;
158  t1->flags[y][x] |= JPEG2000_T1_SIG;
159  if (negative) {
160  t1->flags[y][x + 1] |= JPEG2000_T1_SIG_W | JPEG2000_T1_SGN_W;
161  t1->flags[y][x - 1] |= JPEG2000_T1_SIG_E | JPEG2000_T1_SGN_E;
162  t1->flags[y + 1][x] |= JPEG2000_T1_SIG_N | JPEG2000_T1_SGN_N;
163  t1->flags[y - 1][x] |= JPEG2000_T1_SIG_S | JPEG2000_T1_SGN_S;
164  } else {
165  t1->flags[y][x + 1] |= JPEG2000_T1_SIG_W;
166  t1->flags[y][x - 1] |= JPEG2000_T1_SIG_E;
167  t1->flags[y + 1][x] |= JPEG2000_T1_SIG_N;
168  t1->flags[y - 1][x] |= JPEG2000_T1_SIG_S;
169  }
170  t1->flags[y + 1][x + 1] |= JPEG2000_T1_SIG_NW;
171  t1->flags[y + 1][x - 1] |= JPEG2000_T1_SIG_NE;
172  t1->flags[y - 1][x + 1] |= JPEG2000_T1_SIG_SW;
173  t1->flags[y - 1][x - 1] |= JPEG2000_T1_SIG_SE;
174 }
175 
176 static const uint8_t lut_gain[2][4] = { { 0, 0, 0, 0 }, { 0, 1, 1, 2 } };
177 
179  Jpeg2000CodingStyle *codsty,
180  Jpeg2000QuantStyle *qntsty,
181  int cbps, int dx, int dy,
182  AVCodecContext *avctx)
183 {
184  uint8_t log2_band_prec_width, log2_band_prec_height;
185  int reslevelno, bandno, gbandno = 0, ret, i, j, csize = 1;
186 
187  if (ret=ff_jpeg2000_dwt_init(&comp->dwt, comp->coord, codsty->nreslevels2decode-1, codsty->transform))
188  return ret;
189  for (i = 0; i < 2; i++)
190  csize *= comp->coord[i][1] - comp->coord[i][0];
191 
192  comp->data = av_malloc_array(csize, sizeof(*comp->data));
193  if (!comp->data)
194  return AVERROR(ENOMEM);
195  comp->reslevel = av_malloc_array(codsty->nreslevels, sizeof(*comp->reslevel));
196  if (!comp->reslevel)
197  return AVERROR(ENOMEM);
198  /* LOOP on resolution levels */
199  for (reslevelno = 0; reslevelno < codsty->nreslevels; reslevelno++) {
200  int declvl = codsty->nreslevels - reslevelno; // N_L -r see ISO/IEC 15444-1:2002 B.5
201  Jpeg2000ResLevel *reslevel = comp->reslevel + reslevelno;
202 
203  /* Compute borders for each resolution level.
204  * Computation of trx_0, trx_1, try_0 and try_1.
205  * see ISO/IEC 15444-1:2002 eq. B.5 and B-14 */
206  for (i = 0; i < 2; i++)
207  for (j = 0; j < 2; j++)
208  reslevel->coord[i][j] =
209  ff_jpeg2000_ceildivpow2(comp->coord_o[i][j], declvl - 1);
210  // update precincts size: 2^n value
211  reslevel->log2_prec_width = codsty->log2_prec_widths[reslevelno];
212  reslevel->log2_prec_height = codsty->log2_prec_heights[reslevelno];
213 
214  /* Number of bands for each resolution level */
215  if (reslevelno == 0)
216  reslevel->nbands = 1;
217  else
218  reslevel->nbands = 3;
219 
220  /* Number of precincts wich span the tile for resolution level reslevelno
221  * see B.6 in ISO/IEC 15444-1:2002 eq. B-16
222  * num_precincts_x = |- trx_1 / 2 ^ log2_prec_width) -| - (trx_0 / 2 ^ log2_prec_width)
223  * num_precincts_y = |- try_1 / 2 ^ log2_prec_width) -| - (try_0 / 2 ^ log2_prec_width)
224  * for Dcinema profiles in JPEG 2000
225  * num_precincts_x = |- trx_1 / 2 ^ log2_prec_width) -|
226  * num_precincts_y = |- try_1 / 2 ^ log2_prec_width) -| */
227  if (reslevel->coord[0][1] == reslevel->coord[0][0])
228  reslevel->num_precincts_x = 0;
229  else
230  reslevel->num_precincts_x =
231  ff_jpeg2000_ceildivpow2(reslevel->coord[0][1],
232  reslevel->log2_prec_width) -
233  (reslevel->coord[0][0] >> reslevel->log2_prec_width);
234 
235  if (reslevel->coord[1][1] == reslevel->coord[1][0])
236  reslevel->num_precincts_y = 0;
237  else
238  reslevel->num_precincts_y =
239  ff_jpeg2000_ceildivpow2(reslevel->coord[1][1],
240  reslevel->log2_prec_height) -
241  (reslevel->coord[1][0] >> reslevel->log2_prec_height);
242 
243  reslevel->band = av_malloc_array(reslevel->nbands, sizeof(*reslevel->band));
244  if (!reslevel->band)
245  return AVERROR(ENOMEM);
246 
247  for (bandno = 0; bandno < reslevel->nbands; bandno++, gbandno++) {
248  Jpeg2000Band *band = reslevel->band + bandno;
249  int cblkno, precno;
250  int nb_precincts;
251 
252  /* TODO: Implementation of quantization step not finished,
253  * see ISO/IEC 15444-1:2002 E.1 and A.6.4. */
254  switch (qntsty->quantsty) {
255  uint8_t gain;
256  int numbps;
257  case JPEG2000_QSTY_NONE:
258  /* TODO: to verify. No quantization in this case */
259  band->f_stepsize = 1;
260  break;
261  case JPEG2000_QSTY_SI:
262  /*TODO: Compute formula to implement. */
263  numbps = cbps +
264  lut_gain[codsty->transform == FF_DWT53][bandno + (reslevelno > 0)];
265  band->f_stepsize = SHL(2048 + qntsty->mant[gbandno],
266  2 + numbps - qntsty->expn[gbandno]);
267  break;
268  case JPEG2000_QSTY_SE:
269  /* Exponent quantization step.
270  * Formula:
271  * delta_b = 2 ^ (R_b - expn_b) * (1 + (mant_b / 2 ^ 11))
272  * R_b = R_I + log2 (gain_b )
273  * see ISO/IEC 15444-1:2002 E.1.1 eqn. E-3 and E-4 */
274  /* TODO/WARN: value of log2 (gain_b ) not taken into account
275  * but it works (compared to OpenJPEG). Why?
276  * Further investigation needed. */
277  gain = cbps;
278  band->f_stepsize = pow(2.0, gain - qntsty->expn[gbandno]);
279  band->f_stepsize *= (qntsty->mant[gbandno] / 2048.0 + 1.0);
280  break;
281  default:
282  band->f_stepsize = 0;
283  av_log(avctx, AV_LOG_ERROR, "Unknown quantization format\n");
284  break;
285  }
286  /* FIXME: In openjepg code stespize = stepsize * 0.5. Why?
287  * If not set output of entropic decoder is not correct. */
288  if (!av_codec_is_encoder(avctx->codec))
289  band->f_stepsize *= 0.5;
290 
291  band->i_stepsize = band->f_stepsize * (1 << 16);
292 
293  /* computation of tbx_0, tbx_1, tby_0, tby_1
294  * see ISO/IEC 15444-1:2002 B.5 eq. B-15 and tbl B.1
295  * codeblock width and height is computed for
296  * DCI JPEG 2000 codeblock_width = codeblock_width = 32 = 2 ^ 5 */
297  if (reslevelno == 0) {
298  /* for reslevelno = 0, only one band, x0_b = y0_b = 0 */
299  for (i = 0; i < 2; i++)
300  for (j = 0; j < 2; j++)
301  band->coord[i][j] =
302  ff_jpeg2000_ceildivpow2(comp->coord_o[i][j] - comp->coord_o[i][0],
303  declvl - 1);
304  log2_band_prec_width = reslevel->log2_prec_width;
305  log2_band_prec_height = reslevel->log2_prec_height;
306  /* see ISO/IEC 15444-1:2002 eq. B-17 and eq. B-15 */
307  band->log2_cblk_width = FFMIN(codsty->log2_cblk_width,
308  reslevel->log2_prec_width);
309  band->log2_cblk_height = FFMIN(codsty->log2_cblk_height,
310  reslevel->log2_prec_height);
311  } else {
312  /* 3 bands x0_b = 1 y0_b = 0; x0_b = 0 y0_b = 1; x0_b = y0_b = 1 */
313  /* x0_b and y0_b are computed with ((bandno + 1 >> i) & 1) */
314  for (i = 0; i < 2; i++)
315  for (j = 0; j < 2; j++)
316  /* Formula example for tbx_0 = ceildiv((tcx_0 - 2 ^ (declvl - 1) * x0_b) / declvl) */
317  band->coord[i][j] =
318  ff_jpeg2000_ceildivpow2(comp->coord_o[i][j] - comp->coord_o[i][0] -
319  (((bandno + 1 >> i) & 1) << declvl - 1),
320  declvl);
321  /* TODO: Manage case of 3 band offsets here or
322  * in coding/decoding function? */
323 
324  /* see ISO/IEC 15444-1:2002 eq. B-17 and eq. B-15 */
325  band->log2_cblk_width = FFMIN(codsty->log2_cblk_width,
326  reslevel->log2_prec_width - 1);
327  band->log2_cblk_height = FFMIN(codsty->log2_cblk_height,
328  reslevel->log2_prec_height - 1);
329 
330  log2_band_prec_width = reslevel->log2_prec_width - 1;
331  log2_band_prec_height = reslevel->log2_prec_height - 1;
332  }
333 
334  for (j = 0; j < 2; j++)
335  band->coord[0][j] = ff_jpeg2000_ceildiv(band->coord[0][j], dx);
336  for (j = 0; j < 2; j++)
337  band->coord[1][j] = ff_jpeg2000_ceildiv(band->coord[1][j], dy);
338 
339  band->prec = av_malloc_array(reslevel->num_precincts_x *
340  reslevel->num_precincts_y,
341  sizeof(*band->prec));
342  if (!band->prec)
343  return AVERROR(ENOMEM);
344 
345  nb_precincts = reslevel->num_precincts_x * reslevel->num_precincts_y;
346 
347  for (precno = 0; precno < nb_precincts; precno++) {
348  Jpeg2000Prec *prec = band->prec + precno;
349 
350  /* TODO: Explain formula for JPEG200 DCINEMA. */
351  /* TODO: Verify with previous count of codeblocks per band */
352 
353  /* Compute P_x0 */
354  prec->coord[0][0] = (precno % reslevel->num_precincts_x) *
355  (1 << log2_band_prec_width);
356  prec->coord[0][0] = FFMAX(prec->coord[0][0], band->coord[0][0]);
357 
358  /* Compute P_y0 */
359  prec->coord[1][0] = (precno / reslevel->num_precincts_x) *
360  (1 << log2_band_prec_height);
361  prec->coord[1][0] = FFMAX(prec->coord[1][0], band->coord[1][0]);
362 
363  /* Compute P_x1 */
364  prec->coord[0][1] = prec->coord[0][0] +
365  (1 << log2_band_prec_width);
366  prec->coord[0][1] = FFMIN(prec->coord[0][1], band->coord[0][1]);
367 
368  /* Compute P_y1 */
369  prec->coord[1][1] = prec->coord[1][0] +
370  (1 << log2_band_prec_height);
371  prec->coord[1][1] = FFMIN(prec->coord[1][1], band->coord[1][1]);
372 
373  prec->nb_codeblocks_width =
374  ff_jpeg2000_ceildivpow2(prec->coord[0][1] -
375  prec->coord[0][0],
376  band->log2_cblk_width);
377  prec->nb_codeblocks_height =
378  ff_jpeg2000_ceildivpow2(prec->coord[1][1] -
379  prec->coord[1][0],
380  band->log2_cblk_height);
381 
382  /* Tag trees initialization */
383  prec->cblkincl =
385  prec->nb_codeblocks_height);
386  if (!prec->cblkincl)
387  return AVERROR(ENOMEM);
388 
389  prec->zerobits =
391  prec->nb_codeblocks_height);
392  if (!prec->zerobits)
393  return AVERROR(ENOMEM);
394 
395  prec->cblk = av_malloc_array(prec->nb_codeblocks_width *
396  prec->nb_codeblocks_height,
397  sizeof(*prec->cblk));
398  if (!prec->cblk)
399  return AVERROR(ENOMEM);
400  for (cblkno = 0; cblkno < prec->nb_codeblocks_width * prec->nb_codeblocks_height; cblkno++) {
401  Jpeg2000Cblk *cblk = prec->cblk + cblkno;
402  uint16_t Cx0, Cy0;
403 
404  /* Compute coordinates of codeblocks */
405  /* Compute Cx0*/
406  Cx0 = (prec->coord[0][0] >> band->log2_cblk_width) << band->log2_cblk_width;
407  Cx0 = Cx0 + ((cblkno % prec->nb_codeblocks_width) << band->log2_cblk_width);
408  cblk->coord[0][0] = FFMAX(Cx0, prec->coord[0][0]);
409 
410  /* Compute Cy0*/
411  Cy0 = (prec->coord[1][0] >> band->log2_cblk_height) << band->log2_cblk_height;
412  Cy0 = Cy0 + ((cblkno / prec->nb_codeblocks_width) << band->log2_cblk_height);
413  cblk->coord[1][0] = FFMAX(Cy0, prec->coord[1][0]);
414 
415  /* Compute Cx1 */
416  cblk->coord[0][1] = FFMIN(Cx0 + (1 << band->log2_cblk_width),
417  prec->coord[0][1]);
418 
419  /* Compute Cy1 */
420  cblk->coord[1][1] = FFMIN(Cy0 + (1 << band->log2_cblk_height),
421  prec->coord[1][1]);
422 
423  if((bandno + !!reslevelno) & 1) {
424  cblk->coord[0][0] += comp->reslevel[reslevelno-1].coord[0][1] - comp->reslevel[reslevelno-1].coord[0][0];
425  cblk->coord[0][1] += comp->reslevel[reslevelno-1].coord[0][1] - comp->reslevel[reslevelno-1].coord[0][0];
426  }
427  if((bandno + !!reslevelno) & 2) {
428  cblk->coord[1][0] += comp->reslevel[reslevelno-1].coord[1][1] - comp->reslevel[reslevelno-1].coord[1][0];
429  cblk->coord[1][1] += comp->reslevel[reslevelno-1].coord[1][1] - comp->reslevel[reslevelno-1].coord[1][0];
430  }
431 
432  cblk->zero = 0;
433  cblk->lblock = 3;
434  cblk->length = 0;
435  cblk->lengthinc = 0;
436  cblk->npasses = 0;
437  }
438  }
439  }
440  }
441  return 0;
442 }
443 
445 {
446  int reslevelno, bandno, cblkno, precno;
447  for (reslevelno = 0; reslevelno < codsty->nreslevels; reslevelno++) {
448  Jpeg2000ResLevel *rlevel = comp->reslevel + reslevelno;
449  for (bandno = 0; bandno < rlevel->nbands; bandno++) {
450  Jpeg2000Band *band = rlevel->band + bandno;
451  for(precno = 0; precno < rlevel->num_precincts_x * rlevel->num_precincts_y; precno++) {
452  Jpeg2000Prec *prec = band->prec + precno;
455  for (cblkno = 0; cblkno < prec->nb_codeblocks_width * prec->nb_codeblocks_height; cblkno++) {
456  Jpeg2000Cblk *cblk = prec->cblk + cblkno;
457  cblk->length = 0;
458  cblk->lblock = 3;
459  }
460  }
461  }
462  }
463 }
464 
466 {
467  int reslevelno, bandno, precno;
468  for (reslevelno = 0; reslevelno < codsty->nreslevels; reslevelno++) {
469  Jpeg2000ResLevel *reslevel = comp->reslevel + reslevelno;
470 
471  for (bandno = 0; bandno < reslevel->nbands; bandno++) {
472  Jpeg2000Band *band = reslevel->band + bandno;
473  for (precno = 0; precno < reslevel->num_precincts_x * reslevel->num_precincts_y; precno++) {
474  Jpeg2000Prec *prec = band->prec + precno;
475  av_freep(&prec->zerobits);
476  av_freep(&prec->cblkincl);
477  av_freep(&prec->cblk);
478  }
479 
480  av_freep(&band->prec);
481  }
482  av_freep(&reslevel->band);
483  }
484 
485  ff_dwt_destroy(&comp->dwt);
486  av_freep(&comp->reslevel);
487  av_freep(&comp->data);
488 }