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