FFmpeg: libavcodec/j2k_dwt.c Source File

00001 /*
00002  * Discrete wavelet transform
00003  * Copyright (c) 2007 Kamil Nowosad
00004  *
00005  * This file is part of FFmpeg.
00006  *
00007  * FFmpeg is free software; you can redistribute it and/or
00008  * modify it under the terms of the GNU Lesser General Public
00009  * License as published by the Free Software Foundation; either
00010  * version 2.1 of the License, or (at your option) any later version.
00011  *
00012  * FFmpeg is distributed in the hope that it will be useful,
00013  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00014  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00015  * Lesser General Public License for more details.
00016  *
00017  * You should have received a copy of the GNU Lesser General Public
00018  * License along with FFmpeg; if not, write to the Free Software
00019  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
00020  */
00021 
00028 #include "j2k_dwt.h"
00029 
00030 const static float scale97[] = {1.625786, 1.230174};
00031 
00032 static inline void extend53(int *p, int i0, int i1)
00033 {
00034     p[i0 - 1] = p[i0 + 1];
00035     p[i1    ] = p[i1 - 2];
00036     p[i0 - 2] = p[i0 + 2];
00037     p[i1 + 1] = p[i1 - 3];
00038 }
00039 
00040 static inline void extend97(float *p, int i0, int i1)
00041 {
00042     int i;
00043 
00044     for (i = 1; i <= 4; i++){
00045         p[i0 - i] = p[i0 + i];
00046         p[i1 + i - 1] = p[i1 - i - 1];
00047     }
00048 }
00049 
00050 static void sd_1d53(int *p, int i0, int i1)
00051 {
00052     int i;
00053 
00054     if (i1 == i0 + 1)
00055         return;
00056 
00057     extend53(p, i0, i1);
00058 
00059     for (i = (i0+1)/2 - 1; i < (i1+1)/2; i++)
00060         p[2*i+1] -= (p[2*i] + p[2*i+2]) >> 1;
00061     for (i = (i0+1)/2; i < (i1+1)/2; i++)
00062         p[2*i] += (p[2*i-1] + p[2*i+1] + 2) >> 2;
00063 }
00064 
00065 static void dwt_encode53(DWTContext *s, int *t)
00066 {
00067     int lev,
00068         w = s->linelen[s->ndeclevels-1][0];
00069     int *line = s->linebuf;
00070     line += 3;
00071 
00072     for (lev = s->ndeclevels-1; lev >= 0; lev--){
00073         int lh = s->linelen[lev][0],
00074             lv = s->linelen[lev][1],
00075             mh = s->mod[lev][0],
00076             mv = s->mod[lev][1],
00077             lp;
00078         int *l;
00079 
00080         // HOR_SD
00081         l = line + mh;
00082         for (lp = 0; lp < lv; lp++){
00083             int i, j = 0;
00084 
00085             for (i = 0; i < lh; i++)
00086                 l[i] = t[w*lp + i];
00087 
00088             sd_1d53(line, mh, mh + lh);
00089 
00090             // copy back and deinterleave
00091             for (i =   mh; i < lh; i+=2, j++)
00092                 t[w*lp + j] = l[i];
00093             for (i = 1-mh; i < lh; i+=2, j++)
00094                 t[w*lp + j] = l[i];
00095         }
00096 
00097         // VER_SD
00098         l = line + mv;
00099         for (lp = 0; lp < lh; lp++) {
00100             int i, j = 0;
00101 
00102             for (i = 0; i < lv; i++)
00103                 l[i] = t[w*i + lp];
00104 
00105             sd_1d53(line, mv, mv + lv);
00106 
00107             // copy back and deinterleave
00108             for (i =   mv; i < lv; i+=2, j++)
00109                 t[w*j + lp] = l[i];
00110             for (i = 1-mv; i < lv; i+=2, j++)
00111                 t[w*j + lp] = l[i];
00112         }
00113     }
00114 }
00115 
00116 static void sd_1d97(float *p, int i0, int i1)
00117 {
00118     int i;
00119 
00120     if (i1 == i0 + 1)
00121         return;
00122 
00123     extend97(p, i0, i1);
00124     i0++; i1++;
00125 
00126     for (i = i0/2 - 2; i < i1/2 + 1; i++)
00127         p[2*i+1] -= 1.586134 * (p[2*i] + p[2*i+2]);
00128     for (i = i0/2 - 1; i < i1/2 + 1; i++)
00129         p[2*i] -= 0.052980 * (p[2*i-1] + p[2*i+1]);
00130     for (i = i0/2 - 1; i < i1/2; i++)
00131         p[2*i+1] += 0.882911 * (p[2*i] + p[2*i+2]);
00132     for (i = i0/2; i < i1/2; i++)
00133         p[2*i] += 0.443506 * (p[2*i-1] + p[2*i+1]);
00134 }
00135 
00136 static void dwt_encode97(DWTContext *s, int *t)
00137 {
00138     int lev,
00139         w = s->linelen[s->ndeclevels-1][0];
00140     float *line = s->linebuf;
00141     line += 5;
00142 
00143     for (lev = s->ndeclevels-1; lev >= 0; lev--){
00144         int lh = s->linelen[lev][0],
00145             lv = s->linelen[lev][1],
00146             mh = s->mod[lev][0],
00147             mv = s->mod[lev][1],
00148             lp;
00149         float *l;
00150 
00151         // HOR_SD
00152         l = line + mh;
00153         for (lp = 0; lp < lv; lp++){
00154             int i, j = 0;
00155 
00156             for (i = 0; i < lh; i++)
00157                 l[i] = t[w*lp + i];
00158 
00159             sd_1d97(line, mh, mh + lh);
00160 
00161             // copy back and deinterleave
00162             for (i =   mh; i < lh; i+=2, j++)
00163                 t[w*lp + j] = scale97[mh] * l[i] / 2;
00164             for (i = 1-mh; i < lh; i+=2, j++)
00165                 t[w*lp + j] = scale97[mh] * l[i] / 2;
00166         }
00167 
00168         // VER_SD
00169         l = line + mv;
00170         for (lp = 0; lp < lh; lp++) {
00171             int i, j = 0;
00172 
00173             for (i = 0; i < lv; i++)
00174                 l[i] = t[w*i + lp];
00175 
00176             sd_1d97(line, mv, mv + lv);
00177 
00178             // copy back and deinterleave
00179             for (i =   mv; i < lv; i+=2, j++)
00180                 t[w*j + lp] = scale97[mv] * l[i] / 2;
00181             for (i = 1-mv; i < lv; i+=2, j++)
00182                 t[w*j + lp] = scale97[mv] * l[i] / 2;
00183         }
00184     }
00185 }
00186 
00187 static void sr_1d53(int *p, int i0, int i1)
00188 {
00189     int i;
00190 
00191     if (i1 == i0 + 1)
00192         return;
00193 
00194     extend53(p, i0, i1);
00195 
00196     for (i = i0/2; i < i1/2 + 1; i++)
00197         p[2*i] -= (p[2*i-1] + p[2*i+1] + 2) >> 2;
00198     for (i = i0/2; i < i1/2; i++)
00199         p[2*i+1] += (p[2*i] + p[2*i+2]) >> 1;
00200 }
00201 
00202 static void dwt_decode53(DWTContext *s, int *t)
00203 {
00204     int lev,
00205         w = s->linelen[s->ndeclevels-1][0];
00206     int *line = s->linebuf;
00207     line += 3;
00208 
00209     for (lev = 0; lev < s->ndeclevels; lev++){
00210         int lh = s->linelen[lev][0],
00211             lv = s->linelen[lev][1],
00212             mh = s->mod[lev][0],
00213             mv = s->mod[lev][1],
00214             lp;
00215         int *l;
00216 
00217         // HOR_SD
00218         l = line + mh;
00219         for (lp = 0; lp < lv; lp++){
00220             int i, j = 0;
00221             // copy with interleaving
00222             for (i =   mh; i < lh; i+=2, j++)
00223                 l[i] = t[w*lp + j];
00224             for (i = 1-mh; i < lh; i+=2, j++)
00225                 l[i] = t[w*lp + j];
00226 
00227             sr_1d53(line, mh, mh + lh);
00228 
00229             for (i = 0; i < lh; i++)
00230                 t[w*lp + i] = l[i];
00231         }
00232 
00233         // VER_SD
00234         l = line + mv;
00235         for (lp = 0; lp < lh; lp++){
00236             int i, j = 0;
00237             // copy with interleaving
00238             for (i =   mv; i < lv; i+=2, j++)
00239                 l[i] = t[w*j + lp];
00240             for (i = 1-mv; i < lv; i+=2, j++)
00241                 l[i] = t[w*j + lp];
00242 
00243             sr_1d53(line, mv, mv + lv);
00244 
00245             for (i = 0; i < lv; i++)
00246                 t[w*i + lp] = l[i];
00247         }
00248     }
00249 }
00250 
00251 static void sr_1d97(float *p, int i0, int i1)
00252 {
00253     int i;
00254 
00255     if (i1 == i0 + 1)
00256         return;
00257 
00258     extend97(p, i0, i1);
00259 
00260     for (i = i0/2 - 1; i < i1/2 + 2; i++)
00261         p[2*i] -= 0.443506 * (p[2*i-1] + p[2*i+1]);
00262     for (i = i0/2 - 1; i < i1/2 + 1; i++)
00263         p[2*i+1] -= 0.882911 * (p[2*i] + p[2*i+2]);
00264     for (i = i0/2; i < i1/2 + 1; i++)
00265         p[2*i] += 0.052980 * (p[2*i-1] + p[2*i+1]);
00266     for (i = i0/2; i < i1/2; i++)
00267         p[2*i+1] += 1.586134 * (p[2*i] + p[2*i+2]);
00268 }
00269 
00270 static void dwt_decode97(DWTContext *s, int *t)
00271 {
00272     int lev,
00273         w = s->linelen[s->ndeclevels-1][0];
00274     float *line = s->linebuf;
00275     line += 5;
00276 
00277     for (lev = 0; lev < s->ndeclevels; lev++){
00278         int lh = s->linelen[lev][0],
00279             lv = s->linelen[lev][1],
00280             mh = s->mod[lev][0],
00281             mv = s->mod[lev][1],
00282             lp;
00283         float *l;
00284 
00285         // HOR_SD
00286         l = line + mh;
00287         for (lp = 0; lp < lv; lp++){
00288             int i, j = 0;
00289             // copy with interleaving
00290             for (i =   mh; i < lh; i+=2, j++)
00291                 l[i] = scale97[1-mh] * t[w*lp + j];
00292             for (i = 1-mh; i < lh; i+=2, j++)
00293                 l[i] = scale97[1-mh] * t[w*lp + j];
00294 
00295             sr_1d97(line, mh, mh + lh);
00296 
00297             for (i = 0; i < lh; i++)
00298                 t[w*lp + i] = l[i];
00299         }
00300 
00301         // VER_SD
00302         l = line + mv;
00303         for (lp = 0; lp < lh; lp++){
00304             int i, j = 0;
00305             // copy with interleaving
00306             for (i =   mv; i < lv; i+=2, j++)
00307                 l[i] = scale97[1-mv] * t[w*j + lp];
00308             for (i = 1-mv; i < lv; i+=2, j++)
00309                 l[i] = scale97[1-mv] * t[w*j + lp];
00310 
00311             sr_1d97(line, mv, mv + lv);
00312 
00313             for (i = 0; i < lv; i++)
00314                 t[w*i + lp] = l[i];
00315         }
00316     }
00317 }
00318 
00319 int ff_j2k_dwt_init(DWTContext *s, uint16_t border[2][2], int decomp_levels, int type)
00320 {
00321     int i, j, lev = decomp_levels, maxlen,
00322         b[2][2];
00323 
00324     if ((unsigned)decomp_levels >= FF_DWT_MAX_DECLVLS)
00325         return AVERROR_INVALIDDATA;
00326     s->ndeclevels = decomp_levels;
00327     s->type = type;
00328 
00329     for (i = 0; i < 2; i++)
00330         for(j = 0; j < 2; j++)
00331             b[i][j] = border[i][j];
00332 
00333     maxlen = FFMAX(b[0][1] - b[0][0],
00334                    b[1][1] - b[1][0]);
00335 
00336     while(--lev >= 0){
00337         for (i = 0; i < 2; i++){
00338             s->linelen[lev][i] = b[i][1] - b[i][0];
00339             s->mod[lev][i] = b[i][0] & 1;
00340             for (j = 0; j < 2; j++)
00341                 b[i][j] = (b[i][j] + 1) >> 1;
00342         }
00343     }
00344     if (type == FF_DWT97)
00345         s->linebuf = av_malloc((maxlen + 12) * sizeof(float));
00346     else if (type == FF_DWT53)
00347         s->linebuf = av_malloc((maxlen + 6) * sizeof(int));
00348     else
00349         return -1;
00350 
00351     if (!s->linebuf)
00352         return AVERROR(ENOMEM);
00353 
00354     return 0;
00355 }
00356 
00357 int ff_j2k_dwt_encode(DWTContext *s, int *t)
00358 {
00359     switch(s->type){
00360         case FF_DWT97:
00361             dwt_encode97(s, t); break;
00362         case FF_DWT53:
00363             dwt_encode53(s, t); break;
00364         default:
00365             return -1;
00366     }
00367     return 0;
00368 }
00369 
00370 int ff_j2k_dwt_decode(DWTContext *s, int *t)
00371 {
00372     switch(s->type){
00373         case FF_DWT97:
00374             dwt_decode97(s, t); break;
00375         case FF_DWT53:
00376             dwt_decode53(s, t); break;
00377         default:
00378             return -1;
00379     }
00380     return 0;
00381 }
00382 
00383 void ff_j2k_dwt_destroy(DWTContext *s)
00384 {
00385     av_freep(&s->linebuf);
00386 }