blockDCT2Dpgm.cc

  1 const char *help = "\
  2 progname: blockDCT2Dpgm.cc\n\
  3 code2html: This program reads a pgm image, decomposes it in blocks and applies DCT to each block.\n\
  4 version: Torch3 vision2.0, 2004-2005\n\
  5 (c) Sebastien Marcel (marcel@idiap.ch)\n";
  6 
  7 #include "ImageGray.h"
  8 #include "ipDCT2D.h"
  9 #include "ipEnhanceDCT2D.h"
 10 #include "ipDCT2Dinverse.h"
 11 #include "ipBlock.h"
 12 #include "DiskXFile.h"
 13 #include "CmdLine.h"
 14 
 15 using namespace Torch;
 16 
 17 int main(int argc, char **argv)
 18 {
 19 	char *image_filename;
 20 	int block_size;
 21 	int block_overlap;
 22 	int dimDCT;
 23 	bool verbose;
 24 	bool enhance;
 25 	real lambda;
 26 	bool inverse;
 27 	bool printblocks;
 28   
 29 	
 30   	CmdLine cmd;
 31 	cmd.setBOption("write log", false);
 32   	cmd.info(help);
 33   	cmd.addText("\nArguments:");
 34   	cmd.addSCmdArg("image filename", &image_filename, "image filename");
 35   	cmd.addText("\nOptions:");
 36   	cmd.addICmdOption("-blocksize", &block_size, 8, "size of blocks");
 37   	cmd.addICmdOption("-dim", &dimDCT, 15, "DCT dimensionality");
 38   	cmd.addICmdOption("-blockoverlap", &block_overlap, 4, "overlap between blocks");
 39   	cmd.addBCmdOption("-verbose", &verbose, false, "verbose");
 40   	cmd.addBCmdOption("-enhance", &enhance, false, "enhance");
 41   	cmd.addRCmdOption("-lambda", &lambda, 1.5, "lambda for enhancement");
 42   	cmd.addBCmdOption("-inverse", &inverse, false, "inverse");
 43   	cmd.addBCmdOption("-pb", &printblocks, false, "print blocks");
 44 	cmd.read(argc, argv);
 45 
 46 
 47 	Image *image_in = NULL;
 48 
 49 	image_in = new ImageGray();
 50 	image_in->setBOption("verbose", verbose);
 51 	image_in->load(image_filename);
 52 
 53 	if(verbose)
 54 	{
 55 		print("Image info:\n");
 56 		print("   width = %d\n", image_in->width);
 57 		print("   height = %d\n", image_in->height);
 58 		print("   format = %s (%d)\n", image_in->coding, image_in->n_planes);
 59 	}
 60 
 61 	ipCore *dct2D = NULL;
 62 
 63 	if(enhance) dct2D = new ipEnhanceDCT2D(block_size, "gray", dimDCT, lambda);
 64 	else dct2D = new ipDCT2D(block_size, "gray", dimDCT);
 65 	dct2D->setBOption("verbose", verbose);
 66 
 67 	ipBlock *convolution = NULL;
 68 
 69 	convolution = new ipBlock(image_in->width, image_in->height, "gray", dct2D, dimDCT, block_size, block_overlap);
 70 	convolution->setBOption("verbose", verbose);
 71 	convolution->process(image_in);
 72 
 73 	int convo_output_size = dimDCT;	
 74 	for(int i = 0 ; i < convolution->n_block ; i++)
 75 	{
 76 	   	real *blocks_ = convolution->seq_out->frames[i];
 77 	
 78 	   	if(printblocks)
 79 		{
 80 		   	real energy = 0.0;
 81 			real log_energy;
 82 
 83 			print("%d : [ ", i);
 84 			for(int j = 0 ; j < convo_output_size ; j++)
 85 			{
 86 				print("%g ", blocks_[j]);
 87 				energy += blocks_[j] * blocks_[j];
 88 			}
 89 			energy = sqrt(energy);
 90 			log_energy = energy < 1.0 ? 0.0 : (real) log(energy); 
 91 			print("] (log E=%g)\n", log_energy);
 92 		}
 93 	}
 94 
 95 	if(inverse)
 96 	{
 97 		ipDCT2Dinverse *idct = new ipDCT2Dinverse(dimDCT, block_size);
 98 		
 99 		print("ipDCT (dctdim = %d, output_size=%d) ...\n", dimDCT, convo_output_size);
100 	
101 		real *image_ = new real [image_in->width * image_in->height];
102 		int *contributions_ = new int [image_in->width * image_in->height];
103 
104 		for(int i = 0 ; i < image_in->width * image_in->height ; i++)
105 		{
106 			image_[i] = 0.0;
107 			contributions_[i] = 0;
108 		}
109 
110 
111 		Sequence *seq_tmp = new Sequence(1, convolution->seq_out->frame_size);
112 		
113 		for(int rowblock=0; rowblock<convolution->rowblocks; rowblock++)
114 		{
115 	   		for(int colblock=0; colblock<convolution->colblocks; colblock++) 
116 			{
117 				int row = convolution->row_offset + rowblock*convolution->pixeladv;
118 				int col = convolution->col_offset + colblock*convolution->pixeladv;
119 
120 				// Storing blocks as [colblock][rowblock] (i.e. [x][y])
121 				int index_block = colblock + rowblock * convolution->colblocks;
122 				seq_tmp->copyFrom(convolution->seq_out->frames[index_block]);
123 				
124 				idct->process(seq_tmp);
125    		
126 				real *idct_output = idct->seq_out->frames[0];
127 
128 				int idx_idct_ = 0;
129 				for(int y = row ; y < row + block_size ; y++)
130 				{
131 				   	int bias_ = y * image_in->width;
132 					
133 					for(int x = col ; x < col + block_size ; x++)
134 					{
135 					   	int idx_ = bias_ + x;
136 
137 						image_[idx_] += idct_output[idx_idct_++];
138 						contributions_[idx_]++;
139 					}
140 				}
141 			}
142 		}	
143 
144 		delete seq_tmp;
145 
146 		for(int i = 0 ; i < image_in->width * image_in->height ; i++)
147 			if(contributions_[i] == 0.0) 
148 			{
149 			   	warning("No contributions for that pixel.");
150 
151 			   	image_[i] = 0.0;
152 			}
153 			else 
154 			{
155 			   	if(image_[i] > 255.0)
156 				{
157 			   		image_[i] = 255.0;
158 				}
159 				else if(image_[i] < 0.0)
160 				{
161 			   		image_[i] = 0.0;
162 				}
163 					
164 			   	image_[i] /= (float) contributions_[i];
165 			}
166 			
167 		ImageGray *grayimage_= new ImageGray(image_in->width, image_in->height);
168 	
169 		grayimage_->copyFrom(image_in->width, image_in->height, image_, "gray");
170 		grayimage_->save("inverse.pgm");
171 		
172 		delete grayimage_;
173 		delete [] contributions_;
174 		delete [] image_;
175 		delete idct;
176 	}
177 
178 	delete convolution;
179 	delete dct2D;
180 	delete image_in;
181 
182 	return(0);
183 }