38
41
using namespace vigra;
39
42
using namespace std;
41
47
typedef vigra::BRGBImage::PixelType RGB;
43
void get_gabor_response(string& imagefile, unsigned int& mask, string& model_file, double& threshold,string&
44
mask_format,vector<double>& svm_responses){
46
// Open SVM model file
47
struct svm_model* model;
49
if((model = svm_load_model(model_file.c_str())) == 0){
50
cout << "Couldn't load model file '" << model_file << "'" << endl << endl;
51
for (int j = 0; j < gNumLocs; j++){
52
svm_responses.push_back(0);
56
cout << "Loaded model file:\t" << model_file << endl;
59
// Integers and containers for libsvm
60
int nr_class=svm_get_nr_class(model);
62
struct svm_node *gabor_responses = (struct svm_node *) malloc(max_nr_attr*sizeof(struct svm_node));
63
double *prob_estimates = (double *) malloc(nr_class*sizeof(double));
65
// Open image using Vigra
68
cout << "Opening image file:\t" << imagefile << endl;
70
// Read image given as first argument
71
// File type is determined automatically
72
vigra::ImageImportInfo info(imagefile.c_str());
74
// Create RGB images of appropriate size
75
//vigra::FRGBImage in(info.width(), info.height());
76
vigra::UInt16RGBImage in(info.width(), info.height());
79
importImage(info, destImage(in));
82
double sizefactor = 1;
83
int nw = info.width(), nh = info.height();
85
// In case we want to save filters
86
// Create this directory and change option in Globals.cpp
87
char basename[] = "gabor_filters/celeste";
89
if (info.width() >= info.height()){
90
if (resize_dimension >= info.width() ){
91
resize_dimension = info.width();
94
if (resize_dimension >= info.height()){
95
resize_dimension = info.height();
98
//cout << "Re-size dimenstion:\t" << resize_dimension << endl;
100
cout << "Image dimensions:\t" << info.width() << " x " << info.height() << endl;
102
// Re-size to max dimension
103
if (info.width() > resize_dimension || info.height() > resize_dimension){
105
if (info.width() >= info.height()){
108
sizefactor = (double)resize_dimension/info.width();
110
// calculate new image size
111
nw = resize_dimension;
112
nh = static_cast<int>(0.5 + (sizefactor*info.height()));
115
sizefactor = (double)resize_dimension/info.height();
118
// calculate new image size
119
nw = static_cast<int>(0.5 + (sizefactor*info.width()));
120
nh = resize_dimension;
124
cout << "Scaling by:\t\t" << sizefactor << endl;
125
cout << "New dimensions:\t\t" << nw << " x " << nh << endl;
127
// create a RGB image of appropriate size
128
//vigra::FRGBImage out(nw, nh);
129
vigra::UInt16RGBImage out(nw, nh);
131
// resize the image, using a bi-cubic spline algorithm
132
resizeImageNoInterpolation(srcImageRange(in),destImageRange(out));
137
//exportImage(srcImageRange(in), ImageExportInfo("test.tif").setPixelType("UINT16"));
139
// Convert to LUV colour space
140
UInt16RGBImage luv(in.width(),in.height());
141
transformImage(srcImageRange(in), destImage(luv), RGB2LuvFunctor<double>() );
142
//transformImage(srcImageRange(in), destImage(luv), RGBPrime2LuvFunctor<double>() );
144
//exportImage(srcImageRange(luv), ImageExportInfo("test_luv.tif").setPixelType("UINT16"));
146
// Prepare Gabor API array
147
float *frameBuf = new float[in.width()*in.height()];
148
float *u_values = new float[in.width()*in.height()];
149
float *v_values = new float[in.width()*in.height()];
150
float** pixels = CreateMatrix( (float)0, in.height(), in.width() );
152
// Do something with each pixel...
153
unsigned int counter = 0;
154
vigra::UInt16RGBImage::iterator img_iter(luv.begin()),end(luv.end());
155
//vigra::FRGBImage::iterator img_iter(luv.begin()),end(luv.end());
157
for(; img_iter != end; ++img_iter) {
159
// [0] is L, [1] is U, [2] is V
160
// We only want L for Gabor filtering
161
frameBuf[counter] = (*img_iter)[0];
163
u_values[counter] = (*img_iter)[1];
164
v_values[counter] = (*img_iter)[2];
166
//cout << "Pixel " << counter << " - L: " << (*img_iter)[0] << endl;
167
//cout << "Pixel " << counter << " - U: " << (*img_iter)[1] << endl;
168
//cout << "Pixel " << counter << " - V: " << (*img_iter)[2] << endl;
172
// Prepare framebuf for Gabor API
174
for (int i = 0; i < in.height(); i++ ){
175
for (int j = 0; j < in.width(); j++ ){
176
pixels[i][j] = frameBuf[k];
177
//cout << i << " " << j << " = " << k << " - " << frameBuf[k] << endl;
184
float *response = NULL;
187
// Scale control points by sizefactor
188
for (int j = 0; j < gNumLocs; j++){
190
//cout << sizefactor << ": " << gLocations[j][0] << "," << gLocations[j][1] << " ---> ";
191
gLocations[j][0] = int(gLocations[j][0] * sizefactor);
192
gLocations[j][1] = int(gLocations[j][1] * sizefactor);
193
//cout << gLocations[j][0] << "," << gLocations[j][1] << endl;
195
// Move CPs to border if the filter radius is out of bounds
196
if (gLocations[j][0] <= gRadius){
197
//cout << "Moving CP to border" << endl;
198
gLocations[j][0] = gRadius + 1;
200
if (gLocations[j][1] <= gRadius){
201
//cout << "Moving CP to border" << endl;
202
gLocations[j][1] = gRadius + 1;
204
if (gLocations[j][0] >= nw - gRadius){
205
//cout << "Moving CP to border" << endl;
206
gLocations[j][0] = nw - gRadius - 1;
208
if (gLocations[j][1] >= nh - gRadius){
209
//cout << "Moving CP to border" << endl;
210
gLocations[j][1] = nh - gRadius - 1;
214
// Do Gabor filtering
215
cout << "Generating feature vector by Gabor filtering..." << endl;
216
response = ProcessChannel( pixels, in.height(), in.width(), response, &len, basename );
218
// Turn the response into SVM vector, and add colour features
219
int vector_length = (int)len/gNumLocs;
221
cout << "Classifying control points using SVM..." << endl;
222
for (int j = 0; j < gNumLocs; j++){
224
int pixel_number = gLocations[j][0] + (in.width() * (gLocations[j][1] - 1)) - 1;
225
unsigned int feature = 1;
229
for ( int v = (j * vector_length); v < ((j + 1) * vector_length); v++){
230
gabor_responses[feature-1].index = feature;
231
gabor_responses[feature-1].value = response[v];
232
//cout << feature << ":" << response[v] << " ";
236
// Work out average colour - U + V channels
237
float u_sum = 0, v_sum = 0;
238
unsigned int pixels_in_region = (gRadius * 2)*(gRadius * 2);
240
for (int t = 1 - gRadius; t <= gRadius; t++){
242
unsigned int this_y_pixel = pixel_number + (t * in.width());
244
for (int r = 1 - gRadius; r <= gRadius; r++){
246
unsigned int this_x_pixel = this_y_pixel + r;
247
u_sum += u_values[this_x_pixel];
248
v_sum += v_values[this_x_pixel];
253
float u_ave = (float)u_sum/pixels_in_region;
254
float v_ave = (float)v_sum/pixels_in_region;
256
// Now work out standard deviation for U and V channels
257
u_sum = 0, v_sum = 0;
259
for (int t = 1 - gRadius; t <= gRadius; t++){
261
unsigned int this_y_pixel = pixel_number + (t * in.width());
263
for (int r = 1 - gRadius; r <= gRadius; r++){
265
unsigned int this_x_pixel = this_y_pixel + r;
267
u_sum += pow(u_values[this_x_pixel]-u_ave,2);
268
v_sum += pow(v_values[this_x_pixel]-v_ave,2);
274
float std_u = sqrt(u_sum/(pixels_in_region-1));
275
float std_v = sqrt(v_sum/(pixels_in_region-1));
277
// Add these colour features to feature vector
278
gabor_responses[feature-1].index = feature;
279
gabor_responses[feature-1].value = u_ave;
280
//cout << feature << ":" << u_ave << " ";
282
gabor_responses[feature-1].index = feature;
283
gabor_responses[feature-1].value = std_u;
284
//cout << feature << ":" << std_u << " ";
286
gabor_responses[feature-1].index = feature;
287
gabor_responses[feature-1].value = v_ave;
288
//cout << feature << ":" << v_ave << " ";
290
gabor_responses[feature-1].index = feature;
291
gabor_responses[feature-1].value = std_v;
292
//cout << feature << ":" << std_v << " ";
294
gabor_responses[feature-1].index = feature;
295
gabor_responses[feature-1].value = u_values[pixel_number];
296
//cout << feature << ":" << u_values[pixel_number] << " ";
298
gabor_responses[feature-1].index = feature;
299
gabor_responses[feature-1].value = v_values[pixel_number];
300
//cout << feature << ":" << v_values[pixel_number] << " " << endl;
301
gabor_responses[feature].index = -1;
303
score = svm_predict_probability(model,gabor_responses,prob_estimates);
304
svm_responses.push_back(prob_estimates[0]);
314
// Create mask file name
315
string mask_name = ("");
316
if (imagefile.substr(imagefile.length()-4,1) == (".")){
318
mask_name.append(imagefile.substr(0,imagefile.length()-4));
322
mask_name.append(imagefile.substr(0,imagefile.length()-4));
324
mask_name.append("_mask.");
325
mask_name.append(mask_format);
327
cout << "Generating mask:\t" << mask_name << endl;
328
// Create mask of same dimensions
329
vigra::BRGBImage mask_out(nw, nh);
332
vigra::initImage(srcIterRange(mask_out.upperLeft(),
333
mask_out.upperLeft() + vigra::Diff2D(nw,nh)),
336
float *mask_response = NULL;
340
// Create grid of fiducial points
341
for (int i = gRadius; i < in.height() - gRadius; i += spacing ){
342
for (int j = gRadius; j < in.width() - gRadius; j += spacing ){
345
// Add extra FP at the end of each row in case nw % gRadius
349
// Add extra FP at the end of each row in case nh % gRadius
350
for (int j = gRadius; j < in.width() - gRadius; j += spacing ){
354
// Create the storage matrix
355
gLocations = CreateMatrix( (int)0, gNumLocs, 2);
357
for (int i = gRadius; i < in.height() - gRadius; i += spacing ){
358
for (int j = gRadius; j < in.width() - gRadius; j += spacing ){
360
gLocations[gNumLocs][0] = j;
361
gLocations[gNumLocs][1] = i;
362
//cout << "fPoint " << gNumLocs << ":\t" << i << " " << j << endl;
366
// Add extra FP at the end of each row in case nw % spacing
369
gLocations[gNumLocs][0] = nw - gRadius - 1;
370
gLocations[gNumLocs][1] = i;
371
//cout << "efPoint " << gNumLocs << ":\t" << i << " " << nw - gRadius - 1 << endl;
377
// Add extra FP at the end of each row in case nh % spacing
380
for (int j = gRadius; j < in.width() - gRadius; j += spacing ){
382
gLocations[gNumLocs][0] = j;
383
gLocations[gNumLocs][1] = nh - gRadius - 1;
384
//cout << "efPoint " << gNumLocs << ":\t" << nh - gRadius - 1 << " " << j << endl;
390
//cout << "Total FPs:\t" << gNumLocs << endl;
394
//cout << "Pre-response " << in.height() << ","<< in.width() << " " << gNumLocs << endl;
395
mask_response = ProcessChannel( pixels, in.height(), in.width(), mask_response, &len, basename );
396
//cout << "Post-response" << endl;
398
// Turn the response into SVM vector, and add colour features
400
int vector_length = (int)len/gNumLocs;
402
for ( int j = 0; j < gNumLocs; j++ ){
404
//cout << j << ":"<<endl;
406
unsigned int pixel_number = gLocations[j][0] + (in.width() * (gLocations[j][1] - 1)) - 1;
407
unsigned int feature = 1;
410
// need one more for index = -1
411
if(j >= max_nr_attr - 1){
413
gabor_responses = (struct svm_node *) realloc(gabor_responses,max_nr_attr*sizeof(struct svm_node));
417
for ( int v = (j * vector_length); v < ((j + 1) * vector_length); v++){
418
gabor_responses[feature-1].index = feature;
419
gabor_responses[feature-1].value = mask_response[v];
420
//cout << feature << ":" << mask_response[v] << " ";
424
// Work out average colour - U + V channels
425
float u_sum = 0, v_sum = 0;
426
int pixels_in_region = (gRadius * 2)*(gRadius * 2);
428
for (int t = 1 - gRadius; t <= gRadius; t++){
430
unsigned int this_y_pixel = pixel_number + (t * in.width());
432
for (int r = 1 - gRadius; r <= gRadius; r++){
434
unsigned int this_x_pixel = this_y_pixel + r;
436
u_sum += u_values[this_x_pixel];
437
v_sum += v_values[this_x_pixel];
442
float u_ave = (float)u_sum/pixels_in_region;
443
float v_ave = (float)v_sum/pixels_in_region;
445
// Now work out standard deviation for U and V channels
446
u_sum = 0, v_sum = 0;
448
for (int t = 1 - gRadius; t <= gRadius; t++){
450
unsigned int this_y_pixel = pixel_number + (t * in.width());
452
for (int r = 1 - gRadius; r <= gRadius; r++){
454
unsigned int this_x_pixel = this_y_pixel + r;
456
u_sum += pow(u_values[this_x_pixel]-u_ave,2);
457
v_sum += pow(v_values[this_x_pixel]-v_ave,2);
462
float std_u = sqrt(u_sum/(pixels_in_region-1));
463
float std_v = sqrt(v_sum/(pixels_in_region-1));
465
// Add these colour features to feature vector
466
gabor_responses[feature-1].index = feature;
467
gabor_responses[feature-1].value = u_ave;
468
//cout << feature << ":" << u_ave << " ";
470
gabor_responses[feature-1].index = feature;
471
gabor_responses[feature-1].value = std_u;
472
//cout << feature << ":" << std_u << " ";
474
gabor_responses[feature-1].index = feature;
475
gabor_responses[feature-1].value = v_ave;
476
//cout << feature << ":" << v_ave << " ";
478
gabor_responses[feature-1].index = feature;
479
gabor_responses[feature-1].value = std_v;
480
//cout << feature << ":" << std_v << " ";
482
gabor_responses[feature-1].index = feature;
483
gabor_responses[feature-1].value = u_values[pixel_number];
484
//cout << feature << ":" << u_values[pixel_number] << " ";
486
gabor_responses[feature-1].index = feature;
487
gabor_responses[feature-1].value = v_values[pixel_number];
488
//cout << feature << ":" << v_values[pixel_number] << " " << endl;
489
gabor_responses[feature].index = -1;
491
score = svm_predict_probability(model,gabor_responses,prob_estimates);
492
//cout << score << " " << prob_estimates[0] << endl;
494
if (prob_estimates[0] >= threshold){
496
//cout << "Cloud\t\t(score " << prob_estimates[0] << " > " << threshold << ")" << endl;
498
unsigned int sub_x0 = gLocations[j][0] - gRadius;
499
unsigned int sub_y0 = gLocations[j][1] - gRadius;
500
unsigned int sub_x1 = gLocations[j][0] + gRadius + 1;
501
unsigned int sub_y1 = gLocations[j][1] + gRadius + 1;
503
//cout << sub_x0 << ","<< sub_y0 << " - " << sub_x1 << "," << sub_y1 << endl;
505
// Set region to black
506
vigra::initImage(srcIterRange(mask_out.upperLeft() + vigra::Diff2D(sub_x0, sub_y0),
507
mask_out.upperLeft() + vigra::Diff2D(sub_x1, sub_y1)),
511
//cout << "Non-cloud\t(score " << prob_estimates[0] << " <= " << threshold << ")" << endl;
514
//cout << feature_vector << endl;
518
delete [] mask_response;
520
// Re-size mask to match original image
521
vigra::BRGBImage mask_resize(info.width(),info.height());
522
resizeImageNoInterpolation(srcImageRange(mask_out),destImageRange(mask_resize));
523
exportImage(srcImageRange(mask_resize), ImageExportInfo(mask_name.c_str()).setPixelType("UINT8"));
527
DisposeMatrix( pixels, in.height() );
534
}catch (vigra::StdException & e){
536
// catch any errors that might have occured and print their reason
537
cout << "Unable to open file:\t" << imagefile << endl << endl;
538
cout << e.what() << endl << endl;
539
for (int j = 0; j < gNumLocs; j++){
540
svm_responses.push_back(0);
546
// Free up libsvm stuff
547
free(gabor_responses);
548
free(prob_estimates);
549
svm_destroy_model(model);
50
bool loadSVMmodel(struct svm_model*& model, string& model_file)
52
if((model = svm_load_model(model_file.c_str())) == 0)
54
cout << "Couldn't load model file '" << model_file << "'" << endl << endl;
59
cout << "Loaded model file:\t" << model_file << endl;
65
void destroySVMmodel(struct svm_model*& model)
67
svm_destroy_model(model);
70
// prepare image for use with celeste (downscale, converting in Luv)
71
void prepareCelesteImage(vigra::UInt16RGBImage& rgb,vigra::UInt16RGBImage& luv, int& resize_dimension, double& sizefactor,bool verbose)
80
if (resize_dimension >= nw )
82
resize_dimension = nw;
87
if (resize_dimension >= nh)
89
resize_dimension = nh;
92
//cout << "Re-size dimenstion:\t" << resize_dimension << endl;
94
cout << "Image dimensions:\t" << rgb.width() << " x " << rgb.height() << endl;
96
// Re-size to max dimension
97
vigra::UInt16RGBImage temp;
99
if (rgb.width() > resize_dimension || rgb.height() > resize_dimension)
101
if (rgb.width() >= rgb.height())
103
sizefactor = (double)resize_dimension/rgb.width();
104
// calculate new image size
105
nw = resize_dimension;
106
nh = static_cast<int>(0.5 + (sizefactor*rgb.height()));
110
sizefactor = (double)resize_dimension/rgb.height();
111
// calculate new image size
112
nw = static_cast<int>(0.5 + (sizefactor*rgb.width()));
113
nh = resize_dimension;
118
cout << "Scaling by:\t\t" << sizefactor << endl;
119
cout << "New dimensions:\t\t" << nw << " x " << nh << endl;
122
// create a RGB image of appropriate size
125
// resize the image, using a bi-cubic spline algorithm
126
vigra::resizeImageNoInterpolation(srcImageRange(rgb),destImageRange(temp));
127
//exportImage(srcImageRange(out), ImageExportInfo("test.tif").setPixelType("UINT16"));
132
vigra::copyImage(srcImageRange(rgb),destImage(temp));
135
// Convert to Luv colour space
136
luv.resize(temp.width(),temp.height());
137
vigra::transformImage(srcImageRange(temp), destImage(luv), RGB2LuvFunctor<double>() );
138
//exportImage(srcImageRange(luv), ImageExportInfo("test_luv.tif").setPixelType("UINT16"));
142
// converts the given Luv image into arrays for Gabors filtering
143
void prepareGaborImage(vigra::UInt16RGBImage& luv, float**& pixels)
145
pixels = CreateMatrix( (float)0, luv.height(), luv.width() );
146
// Prepare framebuf for Gabor API, we need only L channel
147
for (int i = 0; i < luv.height(); i++ )
149
for (int j = 0; j < luv.width(); j++ )
151
pixels[i][j] = luv(j,i)[0];
152
//cout << i << " " << j << " = " << k << " - " << frameBuf[k] << endl;
157
//classify the points with SVM
158
vector<double> classifySVM(struct svm_model* model, int gNumLocs,int**& gLocations,int width,int height,int vector_length, float*& response,int gRadius,vigra::UInt16RGBImage& luv,bool needsMoreIndex=false)
160
std::vector<double> svm_response;
161
// Integers and containers for libsvm
162
int max_nr_attr = 56;
163
int nr_class=svm_get_nr_class(model);
164
struct svm_node *gabor_responses = (struct svm_node *) malloc(max_nr_attr*sizeof(struct svm_node));
165
double *prob_estimates = (double *) malloc(nr_class*sizeof(double));
167
for (int j = 0; j < gNumLocs; j++)
169
unsigned int feature = 1;
173
if(j >= max_nr_attr - 1)
176
gabor_responses = (struct svm_node *) realloc(gabor_responses,max_nr_attr*sizeof(struct svm_node));
180
for ( int v = (j * vector_length); v < ((j + 1) * vector_length); v++)
182
gabor_responses[feature-1].index = feature;
183
gabor_responses[feature-1].value = response[v];
184
//cout << feature << ":" << response[v] << " ";
188
// Work out average colour and variance
189
vigra::FindAverageAndVariance<vigra::UInt16RGBImage::PixelType> average;
190
vigra::inspectImage(srcIterRange(
191
luv.upperLeft()+vigra::Diff2D(gLocations[j][0]-gRadius,gLocations[j][1]-gRadius),
192
luv.upperLeft()+vigra::Diff2D(gLocations[j][0]+gRadius,gLocations[j][1]+gRadius)
194
// Add these colour features to feature vector
196
gabor_responses[feature-1].index = feature;
197
gabor_responses[feature-1].value = average.average()[1];
198
//cout << feature << ":" << u_ave << " ";
200
gabor_responses[feature-1].index = feature;
201
gabor_responses[feature-1].value = sqrt(average.variance()[1]);
202
//cout << feature << ":" << std_u << " ";
204
gabor_responses[feature-1].index = feature;
205
gabor_responses[feature-1].value = average.average()[2];
206
//cout << feature << ":" << v_ave << " ";
208
gabor_responses[feature-1].index = feature;
209
gabor_responses[feature-1].value = sqrt(average.variance()[2]);
210
//cout << feature << ":" << std_v << " ";
212
gabor_responses[feature-1].index = feature;
213
gabor_responses[feature-1].value = luv(gLocations[j][0],gLocations[j][1])[1];
214
//cout << feature << ":" << u_values[pixel_number] << " ";
216
gabor_responses[feature-1].index = feature;
217
gabor_responses[feature-1].value = luv(gLocations[j][0],gLocations[j][1])[2];
218
//cout << feature << ":" << v_values[pixel_number] << " " << endl;
219
gabor_responses[feature].index = -1;
221
svm_predict_probability(model,gabor_responses,prob_estimates);
222
svm_response.push_back(prob_estimates[0]);
224
// Free up libsvm stuff
225
free(gabor_responses);
226
free(prob_estimates);
230
// create a grid of control points for creating masks
231
void createGrid(int& gNumLocs,int**& gLocations,int gRadius,int width, int height)
233
int spacing=(gRadius*2)+1;
234
for (int i = gRadius; i < height - gRadius; i += spacing )
236
for (int j = gRadius; j < width - gRadius; j += spacing )
240
// Add extra FP at the end of each row in case width % gRadius
244
// Add extra FP at the end of each row in case nh % gRadius
245
for (int j = gRadius; j < width - gRadius; j += spacing )
250
// Create the storage matrix
251
gLocations = CreateMatrix( (int)0, gNumLocs, 2);
253
for (int i = gRadius; i < height - gRadius; i += spacing )
255
for (int j = gRadius; j < width - gRadius; j += spacing )
257
gLocations[gNumLocs][0] = j;
258
gLocations[gNumLocs][1] = i;
259
//cout << "fPoint " << gNumLocs << ":\t" << i << " " << j << endl;
263
// Add extra FP at the end of each row in case width % spacing
266
gLocations[gNumLocs][0] = width - gRadius - 1;
267
gLocations[gNumLocs][1] = i;
268
//cout << "efPoint " << gNumLocs << ":\t" << i << " " << nw - gRadius - 1 << endl;
273
// Add extra FP at the end of each row in case height % spacing
274
if (height % spacing)
276
for (int j = gRadius; j < width - gRadius; j += spacing )
278
gLocations[gNumLocs][0] = j;
279
gLocations[gNumLocs][1] = height - gRadius - 1;
280
//cout << "efPoint " << gNumLocs << ":\t" << nh - gRadius - 1 << " " << j << endl;
286
//generates the celeste mask on base of given responses and locations
287
void generateMask(vigra::BImage& mask,int& gNumLocs, int**& gLocations,std::vector<double> svm_responses,int gRadius, double threshold)
289
for ( int j = 0; j < gNumLocs; j++ )
291
if (svm_responses[j] >= threshold)
293
unsigned int sub_x0 = gLocations[j][0] - gRadius;
294
unsigned int sub_y0 = gLocations[j][1] - gRadius;
295
unsigned int sub_x1 = gLocations[j][0] + gRadius + 1;
296
unsigned int sub_y1 = gLocations[j][1] + gRadius + 1;
297
//cout << sub_x0 << ","<< sub_y0 << " - " << sub_x1 << "," << sub_y1 << endl;
299
// Set region to black
300
vigra::initImage(srcIterRange(mask.upperLeft() + vigra::Diff2D(sub_x0, sub_y0),
301
mask.upperLeft() + vigra::Diff2D(sub_x1, sub_y1)), 0);
305
//cout << "Non-cloud\t(score " << prob_estimates[0] << " <= " << threshold << ")" << endl;
310
vigra::BImage getCelesteMask(struct svm_model* model, vigra::UInt16RGBImage& input, int radius, float threshold, int resize_dimension,bool adaptThreshold,bool verbose)
312
vigra::UInt16RGBImage luv;
313
double sizefactor=1.0;
314
prepareCelesteImage(input,luv,resize_dimension,sizefactor,verbose);
316
// Prepare Gabor API array
318
prepareGaborImage(luv,pixels);
320
int** gLocations = NULL;
323
// Create grid of fiducial points
324
createGrid(gNumLocs,gLocations,radius,luv.width(),luv.height());
327
float* mask_response=NULL;
328
mask_response = ProcessChannel(pixels,luv.width(),luv.height(),gNumLocs,gLocations,radius,mask_response,&len);
329
// Turn the response into SVM vector, and add colour features
330
vector<double> svm_responses=classifySVM(model,gNumLocs,gLocations,luv.width(),luv.height(),(int)len/gNumLocs,mask_response,radius,luv);
331
delete [] mask_response;
336
for(unsigned int i=0;i<svm_responses.size();i++)
338
if(svm_responses[i]<minVal)
340
minVal=svm_responses[i];
345
threshold=min(minVal+0.1,1.0);
348
// Create mask of same dimensions
349
vigra::BImage mask_out(luv.width(), luv.height(),255);
350
generateMask(mask_out,gNumLocs,gLocations,svm_responses,radius,threshold);
351
// Re-size mask to match original image
352
vigra::BImage mask_resize(input.width(),input.height());
353
resizeImageNoInterpolation(srcImageRange(mask_out),destImageRange(mask_resize));
354
DisposeMatrix(pixels,luv.height());
355
DisposeMatrix(gLocations,gNumLocs);
356
mask_out.resize(0,0);
360
HuginBase::UIntSet getCelesteControlPoints(struct svm_model* model, vigra::UInt16RGBImage& input, HuginBase::CPointVector cps, int radius, float threshold, int resize_dimension,bool verbose)
362
HuginBase::UIntSet cloudCP;
363
vigra::UInt16RGBImage luv;
364
double sizefactor=1.0;
365
prepareCelesteImage(input,luv,resize_dimension,sizefactor,verbose);
367
// Prepare Gabor API array
369
prepareGaborImage(luv,pixels);
371
int gNumLocs = cps.size();
372
int** gLocations = CreateMatrix( (int)0, gNumLocs, 2);
373
for(unsigned int j=0;j<cps.size();j++)
375
HuginBase::ControlPoint cp=cps[j].second;
376
gLocations[j][0] = int(cp.x1 * sizefactor);
377
gLocations[j][1] = int(cp.y1 * sizefactor);
378
// Move CPs to border if the filter radius is out of bounds
379
if (gLocations[j][0] <= radius)
381
gLocations[j][0] = radius + 1;
383
if (gLocations[j][1] <= radius)
385
gLocations[j][1] = radius + 1;
387
if (gLocations[j][0] >= luv.width() - radius)
389
gLocations[j][0] = luv.width() - radius - 1;
391
if (gLocations[j][1] >= luv.height() - radius)
393
gLocations[j][1] = luv.height() - radius - 1;
398
float* response=NULL;
399
response = ProcessChannel(pixels,luv.width(),luv.height(),gNumLocs,gLocations,radius,response,&len);
400
// Turn the response into SVM vector, and add colour features
401
vector<double> svm_responses=classifySVM(model,gNumLocs,gLocations,luv.width(),luv.height(),(int)len/gNumLocs,response,radius,luv);
404
for(unsigned int i=0;i<svm_responses.size();i++)
406
if(svm_responses[i]>=threshold)
408
cloudCP.insert(cps[i].first);
411
DisposeMatrix(pixels,luv.height());
412
DisposeMatrix(gLocations,gNumLocs);
417
} // end of namespace
b'\\ No newline at end of file'
added
removed
src/celeste/Celeste.cpp
