1 /*M///////////////////////////////////////////////////////////////////////////////////////
2 //
3 // IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
4 //
5 // By downloading, copying, installing or using the software you agree to this license.
6 // If you do not agree to this license, do not download, install,
7 // copy or use the software.
8 //
9 //
10 // License Agreement
11 // For Open Source Computer Vision Library
12 //
13 // Copyright (C) 2013, OpenCV Foundation, all rights reserved.
14 // Third party copyrights are property of their respective owners.
15 //
16 // Redistribution and use in source and binary forms, with or without modification,
17 // are permitted provided that the following conditions are met:
18 //
19 // * Redistribution's of source code must retain the above copyright notice,
20 // this list of conditions and the following disclaimer.
21 //
22 // * Redistribution's in binary form must reproduce the above copyright notice,
23 // this list of conditions and the following disclaimer in the documentation
24 // and/or other materials provided with the distribution.
25 //
26 // * The name of the copyright holders may not be used to endorse or promote products
27 // derived from this software without specific prior written permission.
28 //
29 // This software is provided by the copyright holders and contributors "as is" and
30 // any express or implied warranties, including, but not limited to, the implied
31 // warranties of merchantability and fitness for a particular purpose are disclaimed.
32 // In no event shall the Intel Corporation or contributors be liable for any direct,
33 // indirect, incidental, special, exemplary, or consequential damages
34 // (including, but not limited to, procurement of substitute goods or services;
35 // loss of use, data, or profits; or business interruption) however caused
36 // and on any theory of liability, whether in contract, strict liability,
37 // or tort (including negligence or otherwise) arising in any way out of
38 // the use of this software, even if advised of the possibility of such damage.
39 //
40 //M*/
41
42 #include "precomp.hpp"
43 #include "opencv2/photo.hpp"
44 #include <cmath>
45 #include <vector>
46
47 using namespace std;
48 using namespace cv;
49
50 class Decolor
51 {
52 private:
53 Mat kernelx;
54 Mat kernely;
55 int order;
56
57 public:
58 float sigma;
59
60 Decolor();
61 static vector<double> product(const vector <Vec3i> &comb, const double initRGB[3]);
62 double energyCalcu(const vector <double> &Cg, const vector < vector <double> > &polyGrad, const vector <double> &wei) const;
63 void singleChannelGradx(const Mat &img, Mat& dest) const;
64 void singleChannelGrady(const Mat &img, Mat& dest) const;
65 void gradvector(const Mat &img, vector <double> &grad) const;
66 void colorGrad(const Mat &img, vector <double> &Cg) const;
67 static void add_vector(vector <Vec3i> &comb, int &idx, int r,int g,int b);
68 static void add_to_vector_poly(vector < vector <double> > &polyGrad, const vector <double> &curGrad, int &idx1);
69 void weak_order(const Mat &img, vector <double> &alf) const;
70 void grad_system(const Mat &img, vector < vector < double > > &polyGrad,
71 vector < double > &Cg, vector <Vec3i>& comb) const;
72 static void wei_update_matrix(const vector < vector <double> > &poly, const vector <double> &Cg, Mat &X);
73 static void wei_inti(const vector <Vec3i> &comb, vector <double> &wei);
74 void grayImContruct(vector <double> &wei, const Mat &img, Mat &Gray) const;
75 };
76
energyCalcu(const vector<double> & Cg,const vector<vector<double>> & polyGrad,const vector<double> & wei) const77 double Decolor::energyCalcu(const vector <double> &Cg, const vector < vector <double> > &polyGrad, const vector <double> &wei) const
78 {
79 const size_t size = polyGrad[0].size();
80 vector <double> energy(size);
81 vector <double> temp(size);
82 vector <double> temp1(size);
83
84 for(size_t i=0;i< polyGrad[0].size();i++)
85 {
86 double val = 0.0;
87 for(size_t j =0;j<polyGrad.size();j++)
88 val = val + (polyGrad[j][i] * wei[j]);
89 temp[i] = val - Cg[i];
90 temp1[i] = val + Cg[i];
91 }
92
93 for(size_t i=0;i<polyGrad[0].size();i++)
94 energy[i] = -1.0*log(exp(-1.0*pow(temp[i],2)/sigma) + exp(-1.0*pow(temp1[i],2)/sigma));
95
96 double sum = 0.0;
97 for(size_t i=0;i<polyGrad[0].size();i++)
98 sum +=energy[i];
99
100 return (sum/polyGrad[0].size());
101
102 }
103
Decolor()104 Decolor::Decolor()
105 {
106 kernelx = Mat(1,2, CV_32FC1);
107 kernely = Mat(2,1, CV_32FC1);
108 kernelx.at<float>(0,0)=1.0;
109 kernelx.at<float>(0,1)=-1.0;
110 kernely.at<float>(0,0)=1.0;
111 kernely.at<float>(1,0)=-1.0;
112 order = 2;
113 sigma = 0.02f;
114 }
115
product(const vector<Vec3i> & comb,const double initRGB[3])116 vector<double> Decolor::product(const vector <Vec3i> &comb, const double initRGB[3])
117 {
118 vector <double> res(comb.size());
119 for (size_t i=0;i<comb.size();i++)
120 {
121 double dp = 0.0;
122 for(int j=0;j<3;j++)
123 dp = dp + (comb[i][j] * initRGB[j]);
124 res[i] = dp;
125 }
126 return res;
127 }
128
singleChannelGradx(const Mat & img,Mat & dest) const129 void Decolor::singleChannelGradx(const Mat &img, Mat& dest) const
130 {
131 const int w = img.size().width;
132 const Point anchor(kernelx.cols - kernelx.cols/2 - 1, kernelx.rows - kernelx.rows/2 - 1);
133 filter2D(img, dest, -1, kernelx, anchor, 0.0, BORDER_CONSTANT);
134 dest.col(w - 1) = 0.0;
135 }
136
singleChannelGrady(const Mat & img,Mat & dest) const137 void Decolor::singleChannelGrady(const Mat &img, Mat& dest) const
138 {
139 const int h = img.size().height;
140 const Point anchor(kernely.cols - kernely.cols/2 - 1, kernely.rows - kernely.rows/2 - 1);
141 filter2D(img, dest, -1, kernely, anchor, 0.0, BORDER_CONSTANT);
142 dest.row(h - 1) = 0.0;
143 }
144
gradvector(const Mat & img,vector<double> & grad) const145 void Decolor::gradvector(const Mat &img, vector <double> &grad) const
146 {
147 Mat dest;
148 Mat dest1;
149 singleChannelGradx(img,dest);
150 singleChannelGrady(img,dest1);
151
152 Mat d_trans=dest.t();
153 Mat d1_trans=dest1.t();
154
155 const int height = d_trans.size().height;
156 const int width = d_trans.size().width;
157
158 grad.resize(width * height * 2);
159
160 for(int i=0;i<height;i++)
161 for(int j=0;j<width;j++)
162 grad[i*width + j] = d_trans.at<float>(i, j);
163
164 const int offset = width * height;
165 for(int i=0;i<height;i++)
166 for(int j=0;j<width;j++)
167 grad[offset + i * width + j] = d1_trans.at<float>(i, j);
168 }
169
colorGrad(const Mat & img,vector<double> & Cg) const170 void Decolor::colorGrad(const Mat &img, vector <double> &Cg) const
171 {
172 Mat lab;
173
174 cvtColor(img,lab,COLOR_BGR2Lab);
175
176 vector <Mat> lab_channel;
177 split(lab,lab_channel);
178
179 vector <double> ImL;
180 vector <double> Ima;
181 vector <double> Imb;
182
183 gradvector(lab_channel[0],ImL);
184 gradvector(lab_channel[1],Ima);
185 gradvector(lab_channel[2],Imb);
186
187 Cg.resize(ImL.size());
188 for(size_t i=0;i<ImL.size();i++)
189 {
190 const double res = sqrt(pow(ImL[i],2) + pow(Ima[i],2) + pow(Imb[i],2))/100;
191 Cg[i] = res;
192 }
193 }
194
add_vector(vector<Vec3i> & comb,int & idx,int r,int g,int b)195 void Decolor::add_vector(vector <Vec3i> &comb, int &idx, int r,int g,int b)
196 {
197 comb.push_back(Vec3i(r, g, b));
198 idx++;
199 }
200
add_to_vector_poly(vector<vector<double>> & polyGrad,const vector<double> & curGrad,int & idx1)201 void Decolor::add_to_vector_poly(vector < vector <double> > &polyGrad, const vector <double> &curGrad, int &idx1)
202 {
203 polyGrad.push_back(curGrad);
204 idx1++;
205 }
206
weak_order(const Mat & im,vector<double> & alf) const207 void Decolor::weak_order(const Mat &im, vector <double> &alf) const
208 {
209 Mat img;
210 const int h = im.size().height;
211 const int w = im.size().width;
212 if((h + w) > 800)
213 {
214 const double sizefactor = double(800)/(h+w);
215 resize(im, img, Size(cvRound(w*sizefactor), cvRound(h*sizefactor)));
216 }
217 else
218 {
219 img = im;
220 }
221
222 Mat curIm = Mat(img.size(),CV_32FC1);
223 vector <Mat> rgb_channel;
224 split(img,rgb_channel);
225
226 vector <double> Rg, Gg, Bg;
227 gradvector(rgb_channel[2],Rg);
228 gradvector(rgb_channel[1],Gg);
229 gradvector(rgb_channel[0],Bg);
230
231 vector <double> t1(Rg.size()), t2(Rg.size()), t3(Rg.size());
232 vector <double> tmp1(Rg.size()), tmp2(Rg.size()), tmp3(Rg.size());
233
234 const double level = .05;
235
236 for(size_t i=0;i<Rg.size();i++)
237 {
238 t1[i] = (Rg[i] > level) ? 1.0 : 0.0;
239 t2[i] = (Gg[i] > level) ? 1.0 : 0.0;
240 t3[i] = (Bg[i] > level) ? 1.0 : 0.0;
241 tmp1[i] = (Rg[i] < -1.0*level) ? 1.0 : 0.0;
242 tmp2[i] = (Gg[i] < -1.0*level) ? 1.0 : 0.0;
243 tmp3[i] = (Bg[i] < -1.0*level) ? 1.0 : 0.0;
244 }
245
246 alf.resize(Rg.size());
247 for(size_t i =0 ;i < Rg.size();i++)
248 alf[i] = (t1[i] * t2[i] * t3[i]);
249
250 for(size_t i =0 ;i < Rg.size();i++)
251 alf[i] -= tmp1[i] * tmp2[i] * tmp3[i];
252 }
253
grad_system(const Mat & im,vector<vector<double>> & polyGrad,vector<double> & Cg,vector<Vec3i> & comb) const254 void Decolor::grad_system(const Mat &im, vector < vector < double > > &polyGrad,
255 vector < double > &Cg, vector <Vec3i>& comb) const
256 {
257 Mat img;
258 int h = im.size().height;
259 int w = im.size().width;
260 if((h + w) > 800)
261 {
262 const double sizefactor = double(800)/(h+w);
263 resize(im, img, Size(cvRound(w*sizefactor), cvRound(h*sizefactor)));
264 }
265 else
266 {
267 img = im;
268 }
269
270 h = img.size().height;
271 w = img.size().width;
272 colorGrad(img,Cg);
273
274 Mat curIm = Mat(img.size(),CV_32FC1);
275 vector <Mat> rgb_channel;
276 split(img,rgb_channel);
277
278 int idx = 0, idx1 = 0;
279 for(int r=0 ;r <=order; r++)
280 for(int g=0; g<=order;g++)
281 for(int b =0; b <=order;b++)
282 {
283 if((r+g+b)<=order && (r+g+b) > 0)
284 {
285 add_vector(comb,idx,r,g,b);
286 for(int i = 0;i<h;i++)
287 for(int j=0;j<w;j++)
288 curIm.at<float>(i,j)=static_cast<float>(
289 pow(rgb_channel[2].at<float>(i,j),r)*pow(rgb_channel[1].at<float>(i,j),g)*
290 pow(rgb_channel[0].at<float>(i,j),b));
291 vector <double> curGrad;
292 gradvector(curIm,curGrad);
293 add_to_vector_poly(polyGrad,curGrad,idx1);
294 }
295 }
296 }
297
wei_update_matrix(const vector<vector<double>> & poly,const vector<double> & Cg,Mat & X)298 void Decolor::wei_update_matrix(const vector < vector <double> > &poly, const vector <double> &Cg, Mat &X)
299 {
300 const int size = static_cast<int>(poly.size());
301 const int size0 = static_cast<int>(poly[0].size());
302 Mat P = Mat(size, size0, CV_32FC1);
303
304 for (int i = 0; i < size; i++)
305 for (int j = 0; j < size0;j++)
306 P.at<float>(i,j) = static_cast<float>(poly[i][j]);
307
308 const Mat P_trans = P.t();
309 Mat B = Mat(size, size0, CV_32FC1);
310 for(int i =0;i < size;i++)
311 {
312 for(int j = 0, end = int(Cg.size()); j < end;j++)
313 B.at<float>(i,j) = static_cast<float>(poly[i][j] * Cg[j]);
314 }
315
316 Mat A = P*P_trans;
317 solve(A, B, X, DECOMP_NORMAL);
318
319 }
320
wei_inti(const vector<Vec3i> & comb,vector<double> & wei)321 void Decolor::wei_inti(const vector <Vec3i> &comb, vector <double> &wei)
322 {
323 double initRGB[3] = { .33, .33, .33 };
324
325 wei = product(comb,initRGB);
326
327 vector <int> sum(comb.size());
328
329 for(size_t i=0;i<comb.size();i++)
330 sum[i] = (comb[i][0] + comb[i][1] + comb[i][2]);
331
332 for(size_t i=0;i<sum.size();i++)
333 {
334 if(sum[i] == 1)
335 wei[i] = wei[i] * double(1);
336 else
337 wei[i] = wei[i] * double(0);
338 }
339
340 sum.clear();
341
342 }
343
grayImContruct(vector<double> & wei,const Mat & img,Mat & Gray) const344 void Decolor::grayImContruct(vector <double> &wei, const Mat &img, Mat &Gray) const
345 {
346 const int h = img.size().height;
347 const int w = img.size().width;
348
349 vector <Mat> rgb_channel;
350 split(img,rgb_channel);
351
352 int kk =0;
353
354 for(int r =0;r<=order;r++)
355 for(int g=0;g<=order;g++)
356 for(int b=0;b<=order;b++)
357 if((r + g + b) <=order && (r+g+b) > 0)
358 {
359 for(int i = 0;i<h;i++)
360 for(int j=0;j<w;j++)
361 Gray.at<float>(i,j)=Gray.at<float>(i,j) +
362 static_cast<float>(wei[kk])*pow(rgb_channel[2].at<float>(i,j),r)*pow(rgb_channel[1].at<float>(i,j),g)*
363 pow(rgb_channel[0].at<float>(i,j),b);
364
365 kk=kk+1;
366 }
367
368 double minval, maxval;
369 minMaxLoc(Gray, &minval, &maxval);
370
371 Gray -= minval;
372 Gray /= maxval - minval;
373 }
374