libpgf  6.12.24
PGF - Progressive Graphics File
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
WaveletTransform.cpp
Go to the documentation of this file.
00001 /*
00002  * The Progressive Graphics File; http://www.libpgf.org
00003  * 
00004  * $Date: 2006-05-18 16:03:32 +0200 (Do, 18 Mai 2006) $
00005  * $Revision: 194 $
00006  * 
00007  * This file Copyright (C) 2006 xeraina GmbH, Switzerland
00008  * 
00009  * This program is free software; you can redistribute it and/or
00010  * modify it under the terms of the GNU LESSER GENERAL PUBLIC LICENSE
00011  * as published by the Free Software Foundation; either version 2.1
00012  * of the License, or (at your option) any later version.
00013  * 
00014  * This program is distributed in the hope that it will be useful,
00015  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00016  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00017  * GNU General Public License for more details.
00018  * 
00019  * You should have received a copy of the GNU General Public License
00020  * along with this program; if not, write to the Free Software
00021  * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
00022  */
00023 
00028 
00029 #include "WaveletTransform.h"
00030 
00031 #define c1 1    // best value 1
00032 #define c2 2    // best value 2
00033 
00035 // Constructor: Constructs a wavelet transform pyramid of given size and levels.
00036 // @param width The width of the original image (at level 0) in pixels
00037 // @param height The height of the original image (at level 0) in pixels
00038 // @param levels The number of levels (>= 0)
00039 // @param data Input data of subband LL at level 0
00040 CWaveletTransform::CWaveletTransform(UINT32 width, UINT32 height, int levels, DataT* data) 
00041 : m_nLevels(levels + 1)
00042 , m_subband(0) 
00043 {
00044         ASSERT(m_nLevels > 0 && m_nLevels <= MaxLevel + 1);
00045         InitSubbands(width, height, data);
00046 #ifdef __PGFROISUPPORT__
00047         m_ROIindices.SetLevels(levels + 1);
00048 #endif
00049 }
00050 
00052 // Initialize size subbands on all levels
00053 void CWaveletTransform::InitSubbands(UINT32 width, UINT32 height, DataT* data) {
00054         if (m_subband) Destroy();
00055 
00056         // create subbands
00057         m_subband = new CSubband[m_nLevels][NSubbands];
00058 
00059         // init subbands
00060         UINT32 loWidth = width;
00061         UINT32 hiWidth = width;
00062         UINT32 loHeight = height;
00063         UINT32 hiHeight = height;
00064 
00065         for (int level = 0; level < m_nLevels; level++) {
00066                 m_subband[level][LL].Initialize(loWidth, loHeight, level, LL);  // LL
00067                 m_subband[level][HL].Initialize(hiWidth, loHeight, level, HL);  //    HL
00068                 m_subband[level][LH].Initialize(loWidth, hiHeight, level, LH);  // LH
00069                 m_subband[level][HH].Initialize(hiWidth, hiHeight, level, HH);  //    HH
00070                 hiWidth = loWidth >> 1;                 hiHeight = loHeight >> 1;
00071                 loWidth = (loWidth + 1) >> 1;   loHeight = (loHeight + 1) >> 1;
00072         }
00073         if (data) {
00074                 m_subband[0][LL].SetBuffer(data);
00075         }
00076 }
00077 
00079 // Compute fast forward wavelet transform of LL subband at given level and
00080 // stores result on all 4 subbands of level + 1.
00081 // Wavelet transform used in writing a PGF file
00082 // Forward Transform of srcBand and split and store it into subbands on destLevel
00083 // high pass filter at even positions: 1/4(-2, 4, -2)
00084 // low pass filter at odd positions: 1/8(-1, 2, 6, 2, -1)
00085 // @param level A wavelet transform pyramid level (>= 0 && < Levels())
00086 // @param quant A quantization value (linear scalar quantization)
00087 // @return error in case of a memory allocation problem
00088 OSError CWaveletTransform::ForwardTransform(int level, int quant) {
00089         ASSERT(level >= 0 && level < m_nLevels - 1);
00090         const int destLevel = level + 1;
00091         ASSERT(m_subband[destLevel]);
00092         CSubband* srcBand = &m_subband[level][LL]; ASSERT(srcBand);
00093         const UINT32 width = srcBand->GetWidth();
00094         const UINT32 height = srcBand->GetHeight();
00095         DataT* src = srcBand->GetBuffer(); ASSERT(src);
00096         DataT *row0, *row1, *row2, *row3;
00097 
00098         // Allocate memory for next transform level
00099         for (int i=0; i < NSubbands; i++) {
00100                 if (!m_subband[destLevel][i].AllocMemory()) return InsufficientMemory;
00101         }
00102 
00103         if (height >= FilterHeight) {
00104                 // transform LL subband
00105                 // top border handling
00106                 row0 = src; row1 = row0 + width; row2 = row1 + width;
00107                 ForwardRow(row0, width);
00108                 ForwardRow(row1, width);
00109                 ForwardRow(row2, width);
00110                 for (UINT32 k=0; k < width; k++) {
00111                         row1[k] -= ((row0[k] + row2[k] + c1) >> 1);
00112                         row0[k] += ((row1[k] + c1) >> 1);
00113                 }
00114                 LinearToMallat(destLevel, row0, row1, width);
00115                 row0 = row1; row1 = row2; row2 += width; row3 = row2 + width;
00116 
00117                 // middle part
00118                 for (UINT32 i=3; i < height-1; i += 2) {
00119                         ForwardRow(row2, width);
00120                         ForwardRow(row3, width);
00121                         for (UINT32 k=0; k < width; k++) {
00122                                 row2[k] -= ((row1[k] + row3[k] + c1) >> 1);
00123                                 row1[k] += ((row0[k] + row2[k] + c2) >> 2);
00124                         }
00125                         LinearToMallat(destLevel, row1, row2, width);
00126                         row0 = row2; row1 = row3; row2 = row3 + width; row3 = row2 + width;
00127                 }
00128 
00129                 // bottom border handling
00130                 if (height & 1) {
00131                         for (UINT32 k=0; k < width; k++) {
00132                                 row1[k] += ((row0[k] + c1) >> 1);
00133                         }
00134                         LinearToMallat(destLevel, row1, NULL, width);
00135                         row0 = row1; row1 += width;
00136                 } else {
00137                         ForwardRow(row2, width);
00138                         for (UINT32 k=0; k < width; k++) {
00139                                 row2[k] -= row1[k];
00140                                 row1[k] += ((row0[k] + row2[k] + c2) >> 2);
00141                         }
00142                         LinearToMallat(destLevel, row1, row2, width);
00143                         row0 = row1; row1 = row2; row2 += width;
00144                 }
00145         } else {
00146                 // if height is too small
00147                 row0 = src; row1 = row0 + width;
00148                 // first part
00149                 for (UINT32 k=0; k < height; k += 2) {
00150                         ForwardRow(row0, width);
00151                         ForwardRow(row1, width);
00152                         LinearToMallat(destLevel, row0, row1, width);
00153                         row0 += width << 1; row1 += width << 1;
00154                 }
00155                 // bottom
00156                 if (height & 1) {
00157                         LinearToMallat(destLevel, row0, NULL, width);
00158                 }
00159         }
00160 
00161         if (quant > 0) {
00162                 // subband quantization (without LL)
00163                 for (int i=1; i < NSubbands; i++) {
00164                         m_subband[destLevel][i].Quantize(quant);
00165                 }
00166                 // LL subband quantization
00167                 if (destLevel == m_nLevels - 1) {
00168                         m_subband[destLevel][LL].Quantize(quant);
00169                 }
00170         }
00171 
00172         // free source band
00173         srcBand->FreeMemory();
00174         return NoError;
00175 }
00176 
00178 // Forward transform one row
00179 // high pass filter at even positions: 1/4(-2, 4, -2)
00180 // low pass filter at odd positions: 1/8(-1, 2, 6, 2, -1)
00181 void CWaveletTransform::ForwardRow(DataT* src, UINT32 width) {
00182         if (width >= FilterWidth) {
00183                 UINT32 i = 3;
00184 
00185                 // left border handling
00186                 src[1] -= ((src[0] + src[2] + c1) >> 1);
00187                 src[0] += ((src[1] + c1) >> 1);
00188                 
00189                 // middle part
00190                 for (; i < width-1; i += 2) {
00191                         src[i] -= ((src[i-1] + src[i+1] + c1) >> 1);
00192                         src[i-1] += ((src[i-2] + src[i] + c2) >> 2);
00193                 }
00194 
00195                 // right border handling
00196                 if (width & 1) {
00197                         src[i-1] += ((src[i-2] + c1) >> 1);
00198                 } else {
00199                         src[i] -= src[i-1];
00200                         src[i-1] += ((src[i-2] + src[i] + c2) >> 2);
00201                 }
00202         }
00203 }
00204 
00206 // Copy transformed rows loRow and hiRow to subbands LL,HL,LH,HH
00207 void CWaveletTransform::LinearToMallat(int destLevel, DataT* loRow, DataT* hiRow, UINT32 width) {
00208         const UINT32 wquot = width >> 1;
00209         const bool wrem = width & 1;
00210         CSubband &ll = m_subband[destLevel][LL], &hl = m_subband[destLevel][HL];
00211         CSubband &lh = m_subband[destLevel][LH], &hh = m_subband[destLevel][HH];
00212 
00213         if (hiRow) {
00214                 for (UINT32 i=0; i < wquot; i++) {
00215                         ll.WriteBuffer(*loRow++);       // first access, than increment
00216                         hl.WriteBuffer(*loRow++);
00217                         lh.WriteBuffer(*hiRow++);       // first access, than increment
00218                         hh.WriteBuffer(*hiRow++);
00219                 }
00220                 if (wrem) {
00221                         ll.WriteBuffer(*loRow);
00222                         lh.WriteBuffer(*hiRow);
00223                 }
00224         } else {
00225                 for (UINT32 i=0; i < wquot; i++) {
00226                         ll.WriteBuffer(*loRow++);       // first access, than increment
00227                         hl.WriteBuffer(*loRow++);
00228                 }
00229                 if (wrem) ll.WriteBuffer(*loRow);
00230         }
00231 }
00232 
00234 // Compute fast inverse wavelet transform of all 4 subbands of given level and
00235 // stores result in LL subband of level - 1.
00236 // Inverse wavelet transform used in reading a PGF file
00237 // Inverse Transform srcLevel and combine to destBand
00238 // inverse high pass filter for even positions: 1/4(-1, 4, -1)
00239 // inverse low pass filter for odd positions: 1/8(-1, 4, 6, 4, -1)
00240 // @param srcLevel A wavelet transform pyramid level (> 0 && <= Levels())
00241 // @param w [out] A pointer to the returned width of subband LL (in pixels)
00242 // @param h [out] A pointer to the returned height of subband LL (in pixels)
00243 // @param data [out] A pointer to the returned array of image data
00244 // @return error in case of a memory allocation problem
00245 OSError CWaveletTransform::InverseTransform(int srcLevel, UINT32* w, UINT32* h, DataT** data) {
00246         ASSERT(srcLevel > 0 && srcLevel < m_nLevels);
00247         const int destLevel = srcLevel - 1;
00248         ASSERT(m_subband[destLevel]);
00249         CSubband* destBand = &m_subband[destLevel][LL];
00250         UINT32 width, height;
00251 
00252         // allocate memory for the results of the inverse transform 
00253         if (!destBand->AllocMemory()) return InsufficientMemory;
00254         DataT *dest = destBand->GetBuffer(), *origin = dest, *row0, *row1, *row2, *row3;
00255 
00256 #ifdef __PGFROISUPPORT__
00257         PGFRect destROI = destBand->GetROI();   // is valid only after AllocMemory
00258         width = destROI.Width();
00259         height = destROI.Height();
00260         const UINT32 destWidth = width; // destination buffer width
00261         const UINT32 destHeight = height; // destination buffer height
00262 
00263         // update destination ROI
00264         if (destROI.top & 1) {
00265                 destROI.top++;
00266                 origin += destWidth;
00267                 height--;
00268         }
00269         if (destROI.left & 1) {
00270                 destROI.left++;
00271                 origin++;
00272                 width--;
00273         }
00274 
00275         // init source buffer position
00276         for (int i=0; i < NSubbands; i++) {
00277                 UINT32 left = (destROI.left >> 1) - m_subband[srcLevel][i].GetROI().left;
00278                 UINT32 top = (destROI.top >> 1) - m_subband[srcLevel][i].GetROI().top;
00279                 m_subband[srcLevel][i].InitBuffPos(left, top);
00280         }
00281 #else
00282         width = destBand->GetWidth();
00283         height = destBand->GetHeight();
00284         PGFRect destROI(0, 0, width, height);
00285         const UINT32 destWidth = width; // destination buffer width
00286         const UINT32 destHeight = height; // destination buffer height
00287 
00288         // init source buffer position
00289         for (int i=0; i < NSubbands; i++) {
00290                 m_subband[srcLevel][i].InitBuffPos();
00291         }
00292 #endif
00293 
00294         if (destHeight >= FilterHeight) {
00295                 // top border handling
00296                 row0 = origin; row1 = row0 + destWidth;
00297                 MallatToLinear(srcLevel, row0, row1, width);
00298                 for (UINT32 k=0; k < width; k++) {
00299                         row0[k] -= ((row1[k] + c1) >> 1);
00300                 }
00301 
00302                 // middle part
00303                 row2 = row1 + destWidth; row3 = row2 + destWidth;
00304                 for (UINT32 i=destROI.top + 2; i < destROI.bottom - 1; i += 2) {
00305                         MallatToLinear(srcLevel, row2, row3, width);
00306                         for (UINT32 k=0; k < width; k++) {
00307                                 row2[k] -= ((row1[k] + row3[k] + c2) >> 2);
00308                                 row1[k] += ((row0[k] + row2[k] + c1) >> 1);
00309                         }
00310                         InverseRow(row0, width);
00311                         InverseRow(row1, width);
00312                         row0 = row2; row1 = row3; row2 = row1 + destWidth; row3 = row2 + destWidth;
00313                 }
00314 
00315                 // bottom border handling
00316                 if (height & 1) {
00317                         MallatToLinear(srcLevel, row2, NULL, width);
00318                         for (UINT32 k=0; k < width; k++) {
00319                                 row2[k] -= ((row1[k] + c1) >> 1);
00320                                 row1[k] += ((row0[k] + row2[k] + c1) >> 1);
00321                         }
00322                         InverseRow(row0, width);
00323                         InverseRow(row1, width);
00324                         InverseRow(row2, width);
00325                         row0 = row1; row1 = row2; row2 += destWidth;
00326                 } else {
00327                         for (UINT32 k=0; k < width; k++) {
00328                                 row1[k] += row0[k];
00329                         }
00330                         InverseRow(row0, width);
00331                         InverseRow(row1, width);
00332                         row0 = row1; row1 += destWidth;
00333                 }
00334         } else {
00335                 // height is too small
00336                 row0 = origin; row1 = row0 + destWidth;
00337                 // first part
00338                 for (UINT32 k=0; k < height; k += 2) {
00339                         MallatToLinear(srcLevel, row0, row1, width);
00340                         InverseRow(row0, width);
00341                         InverseRow(row1, width);
00342                         row0 += destWidth << 1; row1 += destWidth << 1;
00343                 }
00344                 // bottom
00345                 if (height & 1) {
00346                         MallatToLinear(srcLevel, row0, NULL, width);
00347                         InverseRow(row0, width);
00348                 } 
00349         }
00350 
00351         // free memory of the current srcLevel
00352         for (int i=0; i < NSubbands; i++) {
00353                 m_subband[srcLevel][i].FreeMemory();
00354         }
00355 
00356         // return info
00357         *w = destWidth;
00358         *h = height;
00359         *data = dest;
00360         return NoError;
00361 }
00362 
00364 // Inverse Wavelet Transform of one row
00365 // inverse high pass filter for even positions: 1/4(-1, 4, -1)
00366 // inverse low pass filter for odd positions: 1/8(-1, 4, 6, 4, -1)
00367 void CWaveletTransform::InverseRow(DataT* dest, UINT32 width) {
00368         if (width >= FilterWidth) {
00369                 UINT32 i = 2;
00370 
00371                 // left border handling
00372                 dest[0] -= ((dest[1] + c1) >> 1);
00373 
00374                 // middle part
00375                 for (; i < width - 1; i += 2) {
00376                         dest[i] -= ((dest[i-1] + dest[i+1] + c2) >> 2);
00377                         dest[i-1] += ((dest[i-2] + dest[i] + c1) >> 1);
00378                 }
00379 
00380                 // right border handling
00381                 if (width & 1) {
00382                         dest[i] -= ((dest[i-1] + c1) >> 1);
00383                         dest[i-1] += ((dest[i-2] + dest[i] + c1) >> 1);
00384                 } else {
00385                         dest[i-1] += dest[i-2];
00386                 }
00387         }
00388 }
00389 
00391 // Copy transformed coefficients from subbands LL,HL,LH,HH to interleaved format
00392 void CWaveletTransform::MallatToLinear(int srcLevel, DataT* loRow, DataT* hiRow, UINT32 width) {
00393         const UINT32 wquot = width >> 1;
00394         const bool wrem = width & 1;
00395         CSubband &ll = m_subband[srcLevel][LL], &hl = m_subband[srcLevel][HL];
00396         CSubband &lh = m_subband[srcLevel][LH], &hh = m_subband[srcLevel][HH];
00397 
00398         if (hiRow) {
00399         #ifdef __PGFROISUPPORT__
00400                 const bool storePos = wquot < ll.BufferWidth();
00401                 UINT32 llPos = 0, hlPos = 0, lhPos = 0, hhPos = 0;
00402 
00403                 if (storePos) {
00404                         // save current src buffer positions
00405                         llPos = ll.GetBuffPos(); 
00406                         hlPos = hl.GetBuffPos(); 
00407                         lhPos = lh.GetBuffPos(); 
00408                         hhPos = hh.GetBuffPos(); 
00409                 }
00410         #endif
00411 
00412                 for (UINT32 i=0; i < wquot; i++) {
00413                         *loRow++ = ll.ReadBuffer();// first access, than increment
00414                         *loRow++ = hl.ReadBuffer();// first access, than increment
00415                         *hiRow++ = lh.ReadBuffer();// first access, than increment
00416                         *hiRow++ = hh.ReadBuffer();// first access, than increment
00417                 }
00418 
00419                 if (wrem) {
00420                         *loRow++ = ll.ReadBuffer();// first access, than increment
00421                         *hiRow++ = lh.ReadBuffer();// first access, than increment
00422                 }
00423 
00424         #ifdef __PGFROISUPPORT__
00425                 if (storePos) {
00426                         // increment src buffer positions
00427                         ll.IncBuffRow(llPos); 
00428                         hl.IncBuffRow(hlPos); 
00429                         lh.IncBuffRow(lhPos); 
00430                         hh.IncBuffRow(hhPos); 
00431                 }
00432         #endif
00433 
00434         } else {
00435         #ifdef __PGFROISUPPORT__
00436                 const bool storePos = wquot < ll.BufferWidth();
00437                 UINT32 llPos = 0, hlPos = 0;
00438 
00439                 if (storePos) {
00440                         // save current src buffer positions
00441                         llPos = ll.GetBuffPos(); 
00442                         hlPos = hl.GetBuffPos(); 
00443                 }
00444         #endif
00445 
00446                 for (UINT32 i=0; i < wquot; i++) {
00447                         *loRow++ = ll.ReadBuffer();// first access, than increment
00448                         *loRow++ = hl.ReadBuffer();// first access, than increment
00449                 }
00450                 if (wrem) *loRow++ = ll.ReadBuffer();
00451 
00452         #ifdef __PGFROISUPPORT__
00453                 if (storePos) {
00454                         // increment src buffer positions
00455                         ll.IncBuffRow(llPos); 
00456                         hl.IncBuffRow(hlPos); 
00457                 }
00458         #endif
00459         }
00460 }
00461 
00462 #ifdef __PGFROISUPPORT__
00463 
00464 
00465 
00466 void CWaveletTransform::SetROI(const PGFRect& rect) {
00467         // create tile indices
00468         m_ROIindices.CreateIndices();
00469 
00470         // compute tile indices
00471         m_ROIindices.ComputeIndices(m_subband[0][LL].GetWidth(), m_subband[0][LL].GetHeight(), rect);
00472 
00473         // compute ROIs
00474         UINT32 w, h;
00475         PGFRect r;
00476 
00477         for (int i=0; i < m_nLevels; i++) {
00478                 const PGFRect& indices = m_ROIindices.GetIndices(i);
00479 
00480                 for (int o=0; o < NSubbands; o++) {
00481                         CSubband& subband = m_subband[i][o];
00482 
00483                         subband.SetNTiles(m_ROIindices.GetNofTiles(i)); // must be called before TilePosition()
00484                         subband.TilePosition(indices.left, indices.top, r.left, r.top, w, h);
00485                         subband.TilePosition(indices.right - 1, indices.bottom - 1, r.right, r.bottom, w, h);
00486                         r.right += w;
00487                         r.bottom += h;
00488                         subband.SetROI(r);
00489                 }
00490         }
00491 }
00492 
00494 
00496 void CRoiIndices::CreateIndices() {
00497         if (!m_indices) {
00498                 // create tile indices 
00499                 m_indices = new PGFRect[m_nLevels];
00500         }
00501 }
00502 
00510 void CRoiIndices::ComputeTileIndex(UINT32 width, UINT32 height, UINT32 pos, bool horizontal, bool isMin) {
00511         ASSERT(m_indices);
00512 
00513         UINT32 m;
00514         UINT32 tileIndex = 0;
00515         UINT32 tileMin = 0, tileMax = (horizontal) ? width : height;
00516         ASSERT(pos <= tileMax);
00517 
00518         // compute tile index with binary search
00519         for (int i=m_nLevels - 1; i >= 0; i--) {
00520                 // store values
00521                 if (horizontal) {
00522                         if (isMin) {
00523                                 m_indices[i].left = tileIndex;
00524                         } else {
00525                                 m_indices[i].right = tileIndex + 1;
00526                         }
00527                 } else {
00528                         if (isMin) {
00529                                 m_indices[i].top = tileIndex;
00530                         } else {
00531                                 m_indices[i].bottom = tileIndex + 1;
00532                         }
00533                 }
00534 
00535                 // compute values
00536                 tileIndex <<= 1;
00537                 m = (tileMin + tileMax)/2;
00538                 if (pos >= m) {
00539                         tileMin = m;
00540                         tileIndex++;
00541                 } else {
00542                         tileMax = m;
00543                 }
00544         }
00545 }
00546 
00552 void CRoiIndices::ComputeIndices(UINT32 width, UINT32 height, const PGFRect& rect) {
00553         ComputeTileIndex(width, height, rect.left, true, true);
00554         ComputeTileIndex(width, height, rect.top, false, true);
00555         ComputeTileIndex(width, height, rect.right, true, false);
00556         ComputeTileIndex(width, height, rect.bottom, false, false);
00557 }
00558 
00559 #endif // __PGFROISUPPORT__