PTdecode/CImg-1.3.0/examples/wavelet_atrous.cpp

Tue, 18 Mar 2014 01:27:15 +0000

author
Philip Pemberton <philpem@philpem.me.uk>
date
Tue, 18 Mar 2014 01:27:15 +0000
changeset 23
f2c7acb4a258
parent 5
1204ebf9340d
permissions
-rwxr-xr-x

Update PTdecode to handle output from other Ptouch drivers

     1 /*
     2  #
     3  #  File        : wavelet_atrous.cpp
     4  #                ( C++ source file )
     5  #
     6  #  Description : Performs a 2D or 3D 'a trous' wavelet transform
     7  #                (using a cubic spline) on an image or a video sequence.
     8  #                This file is a part of the CImg Library project.
     9  #                ( http://cimg.sourceforge.net )
    10  #
    11  #  Author      : Renaud Peteri
    12  #                ( Renaud.Peteri(at)mines-paris.org )
    13  #
    14  #  Institution : CWI, Amsterdam
    15  #
    16  #  Date        : February 2005
    17  #
    18  #  References  : Starck, J.-L., Murtagh, F. and Bijaoui, A.,
    19  #                Image Processing and Data Analysis: The Multiscale Approach,
    20  #                Cambridge University Press, 1998.
    21  #                (Hardback and softback, ISBN 0-521-59084-1 and 0-521-59914-8.)
    22  #
    23  #  License     : CeCILL v2.0
    24  #                ( http://www.cecill.info/licences/Licence_CeCILL_V2-en.html )
    25  #
    26  #  This software is governed by the CeCILL  license under French law and
    27  #  abiding by the rules of distribution of free software.  You can  use,
    28  #  modify and/ or redistribute the software under the terms of the CeCILL
    29  #  license as circulated by CEA, CNRS and INRIA at the following URL
    30  #  "http://www.cecill.info".
    31  #
    32  #  As a counterpart to the access to the source code and  rights to copy,
    33  #  modify and redistribute granted by the license, users are provided only
    34  #  with a limited warranty  and the software's author,  the holder of the
    35  #  economic rights,  and the successive licensors  have only  limited
    36  #  liability.
    37  #
    38  #  In this respect, the user's attention is drawn to the risks associated
    39  #  with loading,  using,  modifying and/or developing or reproducing the
    40  #  software by the user in light of its specific status of free software,
    41  #  that may mean  that it is complicated to manipulate,  and  that  also
    42  #  therefore means  that it is reserved for developers  and  experienced
    43  #  professionals having in-depth computer knowledge. Users are therefore
    44  #  encouraged to load and test the software's suitability as regards their
    45  #  requirements in conditions enabling the security of their systems and/or
    46  #  data to be ensured and,  more generally, to use and operate it in the
    47  #  same conditions as regards security.
    48  #
    49  #  The fact that you are presently reading this means that you have had
    50  #  knowledge of the CeCILL license and that you accept its terms.
    51  #
    52 */
    54 #include "CImg.h"
    55 using namespace cimg_library;
    57 // The lines below are necessary when using a non-standard compiler as visualcpp6.
    58 #ifdef cimg_use_visualcpp6
    59 #define std
    60 #endif
    61 #ifdef min
    62 #undef min
    63 #undef max
    64 #endif
    66 #ifndef cimg_imagepath
    67 #define cimg_imagepath "img/"
    68 #endif
    70 CImg<float> mask_x(const unsigned char scale) {
    71   unsigned char d1 = (unsigned char)std::pow(2.0,(double)(scale-1));
    72   unsigned char d2 = (unsigned char)std::pow(2.0,(double)(scale));
    73   unsigned char cx = (unsigned char)std::pow(2.0,(double)(scale));
    74   unsigned char res = (unsigned char)std::pow(2.0,(double)scale);
    75   CImg<float> m(2*res +1,1,1);m.fill(0);
    76   m(cx) = 6.0;
    77   m(cx-d1) =  m(cx+d1) =4.0;
    78   m(cx-d2) =  m(cx+d2) =1.0;
    79   m /= 16.0;
    80   return m;
    81 }
    83 CImg<float> mask_y(const unsigned char scale) {
    84   unsigned char d1 = (unsigned char)std::pow(2.0,(double)(scale-1));
    85   unsigned char d2 = (unsigned char)std::pow(2.0,(double)(scale));
    86   unsigned char cy = (unsigned char)std::pow(2.0,(double)(scale));
    87   unsigned char res = (unsigned char)std::pow(2.0,(double)scale);
    88   CImg<float> m(1,2*res +1);m.fill(0);
    89   m(0,cy) = 6.0;
    90   m(0,cy-d1) =  m(0,cy+d1) =4.0;
    91   m(0,cy-d2) =  m(0,cy+d2) =1.0;
    92   m /= 16.0;
    93   return m;
    94 }
    96 CImg<float> mask_t(const unsigned char scale) {
    97   unsigned char d1 = (unsigned char)std::pow(2.0,(double)(scale-1));
    98   unsigned char d2 = (unsigned char)std::pow(2.0,(double)(scale));
    99   unsigned char ct = (unsigned char)std::pow(2.0,(double)(scale));
   100   unsigned char res = (unsigned char)std::pow(2.0,(double)scale);
   101   CImg<float> m(1,1,2*res +1);m.fill(0);
   102   m(0,0,ct) = 6.0;
   103   m(0,0,ct-d1) =  m(0,0,ct+d1) =4.0;
   104   m(0,0,ct-d2) =  m(0,0,ct+d2) =1.0;
   105   m /= 16.0;
   106   return m;
   107 }
   109 /*------------------
   110   Main procedure
   111   ----------------*/
   112 int main(int argc,char **argv) {
   114   cimg_usage("Perform an 'a trous' wavelet transform (using a cubic spline) on an image or on a video sequence.\n"
   115              "This wavelet transform is undecimated and produces 2 images/videos at each scale. For an example of\n"
   116              "decomposition on a video, try -i img/trees.inr (sequence from the MIT).\n"
   117              "\t(Type -h for help)");
   119   // Read command line parameters
   120   const char
   121     *name_i  = cimg_option("-i",cimg_imagepath "lena.pgm","Input image or video"),
   122     *name_o  = cimg_option("-o","","Name of the multiscale analysis output"),
   123     *axe_dec = cimg_option("-axe",(char*)NULL,"Perform the multiscale decomposition in just one direction ('x', 'y' or 't')");
   124   const unsigned int
   125     s = cimg_option("-s",3,"Scale of decomposition");
   127   const bool help = cimg_option("-h",false,"Display Help");
   128   if(help) exit(0);
   130   // Initialize Image Data
   131   std::fprintf(stderr," - Load image sequence '%s'...\n",cimg::basename(name_i));
   132   const CImg<float> texture_in(name_i);
   133   CImg<float> mask_conv_x, mask_conv_y, mask_conv_t;
   134   CImgList<float> res(s, texture_in.dimx(),texture_in.dimy(),texture_in.dimz());
   135   CImgList<float> wav(s,texture_in.dimx(), texture_in.dimy(), texture_in.dimz());
   136   cimglist_for(res,l) { res(l).fill(0.0); wav(l).fill(0.0);}
   137   unsigned int i;
   139   if (!axe_dec){
   140     // Perform the multiscale decomposition in all directions
   141     for(i=0;i<s;i++){
   142       std::fprintf(stderr," - Performing scale %u ...\n",i+1);
   143       if(i==0){ res(i) =  texture_in;} else {  res(i) = res(i-1);}
   144       mask_conv_x = mask_x(i+1);
   145       res(i) = res(i).get_convolve(mask_conv_x);
   146       mask_conv_y = mask_y(i+1);
   147       res(i) = res(i).get_convolve(mask_conv_y);
   148       mask_conv_t = mask_t(i+1);
   149       res(i) = res(i).get_convolve(mask_conv_t);
   150       if(i==0){wav(i) = texture_in - res(i);}  // res(0) and wav(0) are the 1st scale of decompostion
   151       else {wav(i) = res(i-1) - res(i);}
   152     } }
   154   if (axe_dec) {
   155     // Perform the multiscale decomposition in just one direction
   156     char c;
   157     c = cimg::uncase(axe_dec[0]);
   158     fprintf(stderr," - Decompose the image along axe '%c'\n",c); fflush(stdout);
   160     switch(c) {
   161     case 'x': {
   162       for(i=0;i<s;i++) {
   163         std::fprintf(stderr," - Performing scale %u ...\n",i+1);
   164         if(i==0){ res(i) =  texture_in;} else {  res(i) = res(i-1);}
   165         mask_conv_x = mask_x(i+1);
   166         res(i) = res(i).get_convolve(mask_conv_x);
   167         if(i==0){wav(i) = texture_in - res(i);}
   168         else {wav(i) = res(i-1) - res(i);}}}
   169       break;
   171     case 'y': {
   172       for(i=0;i<s;i++) {
   173         std::fprintf(stderr," - Performing scale %u ...\n",i+1);
   174         if(i==0){ res(i) =  texture_in;} else {  res(i) = res(i-1);}
   175         mask_conv_y = mask_y(i+1);
   176         res(i) = res(i).get_convolve(mask_conv_y);
   177         if(i==0){wav(i) = texture_in - res(i);}
   178         else {wav(i) = res(i-1) - res(i);}}}
   179       break;
   181     case 't': {
   182       for(i=0;i<s;i++) {
   183         std::fprintf(stderr," - Performing scale %u ...\n",i+1);
   184         if(i==0){ res(i) =  texture_in;} else {  res(i) = res(i-1);}
   185         mask_conv_t = mask_t(i+1);
   186         res(i) = res(i).get_convolve(mask_conv_t);
   187         if(i==0){wav(i) = texture_in - res(i);}
   188         else {wav(i) = res(i-1) - res(i);}}}
   189       break;
   191     default: throw CImgException("Error, unknow decompostion axe '%c', try 'x', 'y' or 't'",c);
   192     }
   193     fputc('\n',stderr);
   194   }
   196   if (*name_o){
   197     // Save the Multi-Scale Analysis
   198     std::fprintf(stderr," - Saving of all output sequences : %s in the msa/ directory... \n",cimg::basename(name_o));
   199     int count = 1; // res0 = original image
   200     char filename[256] = "", filename_wav[256] = "";
   201     char STmp[3] = "";
   202     system("mkdir msa");
   203     for(i=0;i<s;i++){
   204       strcpy( filename, "msa/res" );
   205       strcpy( filename_wav, "msa/wav" );
   206       if( count < 10 )
   207         { strcat( filename, "0" );strcat( filename_wav, "0" );}
   208       sprintf( STmp, "%d_", count );
   209       strcat( filename, STmp ); strcat( filename_wav, STmp );
   210       strcat( filename,name_o);strcat( filename_wav,name_o);
   211       res(i).save(filename);
   212       wav(i).save(filename_wav);
   213       count++;
   214     }
   215   }
   217   // Result visualization
   218   const float value = 255;
   219   for(i=0;i<s;i++) {
   220     res[i].normalize(0,255).draw_text(2,2,"Scale %d",&value,0,1,11,i);
   221     wav[i].normalize(0,255).draw_text(2,2,"Scale %d",&value,0,1,11,i);
   222   }
   224   CImgDisplay disp(res,"Approximations levels by increasing scale",0);
   225   CImgDisplay disp2(wav,"Wavelet coefficients by increasing scale",0);
   226   while ( !disp.is_closed  && !disp.is_keyQ  && !disp.is_keyESC &&
   227           !disp2.is_closed && !disp2.is_keyQ && !disp2.is_keyESC ) {
   228     if (disp.is_resized) disp.resize().display(res);
   229     if (disp2.is_resized) disp2.resize().display(wav);
   230     CImgDisplay::wait(disp,disp2);
   231   }
   233   return 0;
   234 }