/*************************************************************************** * Copyright (C) 2008 by Tim Nugent * * timnugent@gmail.com * * * * This program is free software; you can redistribute it and/or modify * * it under the terms of the GNU General Public License as published by * * the Free Software Foundation; either version 2 of the License, or * * (at your option) any later version. * * * * This program is distributed in the hope that it will be useful, * * but WITHOUT ANY WARRANTY; without even the implied warranty of * * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * * GNU General Public License for more details. * * * * You should have received a copy of the GNU General Public License * * along with this program; if not, write to the * * Free Software Foundation, Inc., * * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. * ***************************************************************************/ #include #include "vigra/stdimage.hxx" #include "vigra/resizeimage.hxx" #include "vigra/impex.hxx" #include "vigra/colorconversions.hxx" #include #include #include #include "Gabor.h" #include "Utilities.h" #include "CelesteGlobals.h" using namespace vigra; using namespace std; void usage (char* filename){ cout << "Usage: " << filename << " " << endl; cout << "or for a positive example: " << filename << " -p " << endl; cout << "or for a negative example: " << filename << " -n " << endl; cout << "(supported formats: " << vigra::impexListFormats() << ")" << endl; } int main(int argc, char ** argv){ unsigned int a = 1; unsigned int pos = 0; unsigned int neg = 0; string image_file,mask_file; vigra::FRGBImage mask_in; vigra::FRGBImage::iterator mask_iter; // Deal with arguments if(argc < 3 || argc >= 4){ usage(argv[0]); exit(1); } while( a < argc ){ if( argv[a][0] == '-'){ a++; switch(argv[a-1][1]){ // Put some more argument options in here later case 'h' : {usage(argv[0]); exit(1);} case 'p' : {pos = 1; cerr << "Ignoring mask. Image is a positive example." << endl; image_file = argv[argc-1]; break; } case 'n' : {neg = 1; cerr << "Ignoring mask. Image is a negative example." << endl; image_file = argv[argc-1]; break; } default : {usage(argv[0]); exit(1);} } } a++; } if (pos == 0 && neg == 0){ image_file = argv[argc-2]; mask_file = argv[argc-1]; //cout << "Image is " << argv[argc-2] << endl; //cout << "Mask is " << argv[argc-1] << endl; } try{ cerr << "Opening image file: " << image_file << endl; // Read image given as first argument // File type is determined automatically vigra::ImageImportInfo info(image_file.c_str()); // Create RGB images of appropriate size vigra::FRGBImage in(info.width(), info.height()); // Import the image importImage(info, destImage(in)); // Deal with mask unless -p || -n flags are on if (pos == 0 && neg == 0){ cerr << "Opening mask file: " << mask_file << endl; // Read mask given as second argument vigra::ImageImportInfo mask_info(mask_file.c_str()); // Check dimensions are the same if (info.width() != mask_info.width() || info.height() != mask_info.height()){ cerr << "Error - image and mask files are different dimension:" << endl; cerr << "Image:\t" << info.width() << " x " << info.height() << endl; cerr << "Mask:\t" << mask_info.width() << " x " << mask_info.height() << endl << endl; return 1; } // Create RGB images of appropriate size vigra::FRGBImage tmp_mask_in(info.width(), info.height()); // Import the mask importImage(mask_info, destImage(tmp_mask_in)); mask_in = tmp_mask_in; } // Max dimension const int resize_dimension = 800; // Re-size to 800 max dimension if (info.width() >= resize_dimension || info.height() >= resize_dimension){ // Method (0 - pixel repetition, 1 - linear, 2 - spline // Pixel repetition seems to be the fastest int method = 0; double sizefactor; if (info.width() >= info.height()){ sizefactor = (double)resize_dimension/info.width(); }else{ sizefactor = (double)resize_dimension/info.height(); } // calculate new image size int nw = (int)(sizefactor*info.width()); int nh = (int)(sizefactor*info.height()); cerr << "Image dimensions:\t" << info.width() << " x " << info.height() << endl; cerr << "Scaling by:\t\t" << sizefactor << endl; cerr << "New dimensions:\t\t" << nw << " x " << nh << endl; // create a RGB image of appropriate size vigra::FRGBImage out(nw, nh); switch(method){ case 0: // resize the image, using a bi-cubic spline algorithms resizeImageNoInterpolation(srcImageRange(in),destImageRange(out)); break; case 1: // resize the image, using a bi-cubic spline algorithms resizeImageLinearInterpolation(srcImageRange(in),destImageRange(out)); break; default: // resize the image, using a bi-cubic spline algorithms resizeImageSplineInterpolation(srcImageRange(in),destImageRange(out)); } in = out; // Deal with mask unless -p || -n flags are on if (pos == 0 && neg == 0){ vigra::FRGBImage mask_out(nw, nh); switch(method){ case 0: // resize the image, using a bi-cubic spline algorithms resizeImageNoInterpolation(srcImageRange(mask_in),destImageRange(mask_out)); break; case 1: // resize the image, using a bi-cubic spline algorithms resizeImageLinearInterpolation(srcImageRange(mask_in),destImageRange(mask_out)); break; default: // resize the image, using a bi-cubic spline algorithms resizeImageSplineInterpolation(srcImageRange(mask_in),destImageRange(mask_out)); } mask_in = mask_out; } }else{ cerr << "Image dimensions:\t" << info.width() << " x " << info.height() << endl; } cerr << "Total pixels:\t\t" << in.width()*in.height() << endl; // Convert to LUV colour space cerr << "Converting to LUV colourspace..." << endl; FRGBImage luv(in.width(),in.height()); transformImage(srcImageRange(in), destImage(luv), RGBPrime2LuvFunctor() ); // Create grid of fiducial points cerr << "Generating fiducial points..." << endl; for (int i = gRadius; i < in.height() - gRadius; i += spacing ){ for (int j = gRadius; j < in.width() - gRadius; j += spacing ){ gNumLocs++; } } // Create the storage matrix gLocations = CreateMatrix( (int)0, gNumLocs, 2); gNumLocs = 0; for (int i = gRadius; i < in.height() - gRadius; i += spacing ){ for (int j = gRadius; j < in.width() - gRadius; j += spacing ){ gLocations[gNumLocs][0] = j; gLocations[gNumLocs][1] = i; //cerr << "fPoint " << gNumLocs << ":\t" << i << " " << j << endl; gNumLocs++; } } cerr << "Total fiducial points: " << gNumLocs << endl; // Prepare Gabor API array float *frameBuf = new float[in.width()*in.height()]; float *u_values = new float[in.width()*in.height()]; float *v_values = new float[in.width()*in.height()]; float** pixels = CreateMatrix( (float)0, in.height(), in.width() ); // Do something with each pixel... int counter = 0; vigra::FRGBImage::iterator img_iter(luv.begin()),end(luv.end()); for(; img_iter != end; ++img_iter) { // [0] is L, [1] is U, [2] is V // We only want L for Gabor filtering frameBuf[counter] = (*img_iter)[0]; u_values[counter] = (*img_iter)[1]; v_values[counter] = (*img_iter)[2]; //cout << "Pixel " << counter << " - L: " << (*img_iter)[0] << endl; //cout << "Pixel " << counter << " - U: " << (*img_iter)[1] << endl; //cout << "Pixel " << counter << " - V: " << (*img_iter)[2] << endl; counter++; } int k = 0; for (int i = 0; i < in.height(); i++ ){ for (int j = 0; j < in.width(); j++ ){ pixels[i][j] = frameBuf[k]; //cout << i << " " << j << " = " << k << " - " << frameBuf[k] << endl; k++; } } char basename[] = "gabor_filters/celeste"; float *response = NULL; int len = 0; response = ProcessChannel( pixels, in.height(), in.width(), response, &len, basename ); // Deal with mask unless -p || -n flags are on if (pos == 0 && neg == 0){ // Iterator for mask vigra::FRGBImage::iterator tmp_mask_iter(mask_in.begin()),mask_end(mask_in.end()); mask_iter = tmp_mask_iter; mask_end = NULL; } cerr << "Generating feature vector by Gabor filtering..." << endl; // Do something nice with response int vector_length = (int)len/gNumLocs; for ( int j = 0; j < gNumLocs; j++ ){ //cerr << gLocations[j][0] << "," << gLocations[j][1] << endl; int pixel_number = gLocations[j][0] + (in.width() * (gLocations[j][1] - 1)) - 1; // Apply label if (pos){ cout << "+1 "; }else if (neg){ cout << "-1 "; }else{ // In the mask, clouds are black, non-clouds are white if (mask_iter[pixel_number][0] != 0 || mask_iter[pixel_number][1] != 0 || mask_iter[pixel_number][2] != 0){ // So this must be non-cloud = -ve example cout << "-1 "; }else{ // So this must be cloud = +ve example cout << "+1 "; } } int feature = 1; for ( int v = (j * vector_length); v < ((j + 1) * vector_length); v++){ cout << feature << ":" << response[v] << " "; feature++; } // Work out average colour - U + V channels float u_sum = 0, v_sum = 0; int pixels_in_region = (gRadius * 2)*(gRadius * 2); // height //cout << endl << "Pixel " << pixel_number << " x:" << gLocations[j][0] << " y:" << //gLocations[j][1] << " width: " << in.width() << endl; for (int t = 1 - gRadius; t <= gRadius; t++){ int this_y_pixel = pixel_number + (t * in.width()); //cout << "t = " << t << " pixel= " << pixel_number << " + " << (t * in.width()) << endl; // width for (int r = 1 - gRadius; r <= gRadius; r++){ int this_x_pixel = this_y_pixel + r; //if (this_x_pixel == pixel_number){ // cout << "*" << this_x_pixel << "* "; //}else{ // cout << this_x_pixel << " "; //} u_sum += u_values[this_x_pixel]; v_sum += v_values[this_x_pixel]; //cout << "u = " << u_values[this_x_pixel] << "(" << u_sum << ")" << endl; //cout << "v = " << v_values[this_x_pixel] << "(" << v_sum << ")" << endl; } //cout << endl; } //cout << endl; float u_ave = (float)u_sum/pixels_in_region; float v_ave = (float)v_sum/pixels_in_region; // Now work out standard deviation for U and V channels u_sum = 0, v_sum = 0; for (int t = 1 - gRadius; t <= gRadius; t++){ int this_y_pixel = pixel_number + (t * in.width()); for (int r = 1 - gRadius; r <= gRadius; r++){ int this_x_pixel = this_y_pixel + r; u_sum += pow(u_values[this_x_pixel]-u_ave,2); v_sum += pow(v_values[this_x_pixel]-v_ave,2); } } float std_u = sqrt(u_sum/(pixels_in_region-1)); float std_v = sqrt(v_sum/(pixels_in_region-1)); //cerr << "U sum = " << u_sum << " U ave = " << u_ave << " U std = " << std_u << endl; //cerr << "V sum = " << v_sum << " V ave = " << v_ave << " V std = " << std_v << endl << endl; cout << feature << ":" << u_ave << " "; feature++; cout << feature << ":" << std_u << " "; feature++; cout << feature << ":" << v_ave << " "; feature++; cout << feature << ":" << std_v << " "; feature++; // Add the U/V values for this actual pixel cout << feature << ":" << u_values[pixel_number] << " "; feature++; cout << feature << ":" << v_values[pixel_number] << " "; feature++; //cout << "# FP_" << j << "_Pixel_" << pixel_number << "_Pos_" << gLocations[j][0] << "-" << gLocations[j][1] << //"_Mask_RGB_" << mask_iter[pixel_number][0] << "," << mask_iter[pixel_number][1] << "," << //mask_iter[pixel_number][2] << endl; cout << endl; } delete[] response; delete[] frameBuf; delete[] u_values; delete[] v_values; // Export images //exportImage(srcImageRange(in), ImageExportInfo("output_images/image_reduced.jpg")); //exportImage(srcImageRange(luv), ImageExportInfo("output_images/luv_reduced.jpg")); // //if (pos == 0 && neg == 0){ // exportImage(srcImageRange(mask_in), ImageExportInfo("output_images/mask_reduced.jpg")); //} //cout << "Done." << endl; }catch (vigra::StdException & e){ // catch any errors that might have occured and print their reason std::cout << e.what() << std::endl; return 1; } return 0; }