API tools faq
paste
Advertisement
Guest User

Untitled

a guest
Feb 19th, 2019
74
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 28.35 KB | None | 0 0
raw download clone embed print report
  1. #define Xxx ixxG
  2.  
  3. Tools::Tools()
  4. {
  5.  
  6. }
  7.  
  8.  
  9. float gF(float val){
  10. return 1/(1+(val*val)/16);
  11. }
  12. float fF(float val){
  13. return exp(-val);
  14. }
  15. float NF(float x,float y){
  16. return sqrt(x*x+y*y);
  17. }
  18.  
  19. Mat Tools::_Convolution(Mat image,Mat kernel){
  20. Mat dst;
  21. image.convertTo(image,CV_64F,1.0/255.0);
  22. filter2D(image, dst, -1 , kernel, Point( -1, -1 ), 0, BORDER_DEFAULT );
  23. image.convertTo(image,CV_8UC1,255.0);
  24. return dst;
  25. }
  26.  
  27. Mat Tools::_DiZenzoGradient(Mat image){
  28. Mat dst;
  29. Mat dstChannel[3];
  30. int progressMCM;
  31.  
  32. image.copyTo(dst);
  33.  
  34. double ix,iy;
  35.  
  36. Mat ixB,iyB;
  37. Mat ixG,iyG;
  38. Mat ixR,iyR;
  39.  
  40. double g11,g12,g22;
  41. double lambdaPlus,lambdaMin;
  42.  
  43. split(dst,dstChannel);
  44.  
  45.  
  46. ixB = Tools::_Convolution(dstChannel[0],Tools::_Ix());
  47. iyB = Tools::_Convolution(dstChannel[0],Tools::_Iy());
  48.  
  49. ixG = Tools::_Convolution(dstChannel[1],Tools::_Ix());
  50. iyG = Tools::_Convolution(dstChannel[1],Tools::_Iy());
  51.  
  52. ixR = Tools::_Convolution(dstChannel[2],Tools::_Ix());
  53. iyR = Tools::_Convolution(dstChannel[2],Tools::_Iy());
  54.  
  55.  
  56. Mat cM = Mat::zeros(Size(image.cols, image.rows), CV_64F);
  57.  
  58. for(int i = 1; i < image.rows - 1; i++){
  59. for(int j = 1; j < image.cols - 1; j++){
  60.  
  61. g11= pow(ixR.at<double>(i,j),2)+pow(ixG.at<double>(i,j),2)+pow(ixB.at<double>(i,j),2);
  62. g12= ixR.at<double>(i,j)*iyR.at<double>(i,j)+ixG.at<double>(i,j)*iyG.at<double>(i,j)+ixB.at<double>(i,j)*iyB.at<double>(i,j);
  63. g22= pow(iyR.at<double>(i,j),2)+pow(iyG.at<double>(i,j),2)+pow(iyB.at<double>(i,j),2);
  64.  
  65. lambdaPlus=( g11 + g22 + sqrt( pow( (g11-g22),2 )+4*pow( g12,2 ) ) )/2;
  66. lambdaMin=( g11 + g22 - sqrt( pow( (g11-g22),2 )+4*pow( g12,2 ) ) )/2;
  67.  
  68. cM.at<double>(i,j) = sqrt(lambdaPlus);
  69.  
  70. }
  71. }
  72.  
  73. cM.convertTo(cM, CV_8UC1, 255.0/1.0);
  74. return cM;
  75. }
  76.  
  77. Mat Tools::_SapiroGradient(Mat image){
  78. Mat dst;
  79. Mat dstChannel[3];
  80. int progressMCM;
  81.  
  82. image.copyTo(dst);
  83.  
  84. double ix,iy;
  85.  
  86. Mat ixB,iyB;
  87. Mat ixG,iyG;
  88. Mat ixR,iyR;
  89.  
  90. double g11,g12,g22;
  91. double lambdaPlus,lambdaMin;
  92.  
  93. split(dst,dstChannel);
  94.  
  95.  
  96. ixB = Tools::_Convolution(dstChannel[0],Tools::_Ix());
  97. iyB = Tools::_Convolution(dstChannel[0],Tools::_Iy());
  98.  
  99. ixG = Tools::_Convolution(dstChannel[1],Tools::_Ix());
  100. iyG = Tools::_Convolution(dstChannel[1],Tools::_Iy());
  101.  
  102. ixR = Tools::_Convolution(dstChannel[2],Tools::_Ix());
  103. iyR = Tools::_Convolution(dstChannel[2],Tools::_Iy());
  104.  
  105.  
  106. Mat cM = Mat::zeros(Size(image.cols, image.rows), CV_64F);
  107.  
  108. for(int i = 1; i < image.rows - 1; i++){
  109. for(int j = 1; j < image.cols - 1; j++){
  110.  
  111. g11= pow(ixR.at<double>(i,j),2)+pow(ixG.at<double>(i,j),2)+pow(ixB.at<double>(i,j),2);
  112. g12= ixR.at<double>(i,j)*iyR.at<double>(i,j)+ixG.at<double>(i,j)*iyG.at<double>(i,j)+ixB.at<double>(i,j)*iyB.at<double>(i,j);
  113. g22= pow(iyR.at<double>(i,j),2)+pow(iyG.at<double>(i,j),2)+pow(iyB.at<double>(i,j),2);
  114.  
  115. lambdaPlus=( g11 + g22 + sqrt( pow( (g11-g22),2 )+4*pow( g12,2 ) ) )/2;
  116. lambdaMin=( g11 + g22 - sqrt( pow( (g11-g22),2 )+4*pow( g12,2 ) ) )/2;
  117.  
  118. cM.at<double>(i,j) = sqrt(lambdaPlus-lambdaMin);
  119.  
  120. }
  121. }
  122.  
  123. cM.convertTo(cM, CV_8UC1, 255.0/1.0);
  124. return cM;
  125. }
  126.  
  127. Mat Tools::_BlomgrenGradient(Mat image){
  128. Mat dst;
  129. Mat dstChannel[3];
  130. int progressMCM;
  131.  
  132. image.copyTo(dst);
  133.  
  134. double ix,iy;
  135.  
  136. Mat ixB,iyB;
  137. Mat ixG,iyG;
  138. Mat ixR,iyR;
  139.  
  140. double g11,g12,g22;
  141. double lambdaPlus,lambdaMin;
  142.  
  143. split(dst,dstChannel);
  144.  
  145.  
  146. ixB = Tools::_Convolution(dstChannel[0],Tools::_Ix());
  147. iyB = Tools::_Convolution(dstChannel[0],Tools::_Iy());
  148.  
  149. ixG = Tools::_Convolution(dstChannel[1],Tools::_Ix());
  150. iyG = Tools::_Convolution(dstChannel[1],Tools::_Iy());
  151.  
  152. ixR = Tools::_Convolution(dstChannel[2],Tools::_Ix());
  153. iyR = Tools::_Convolution(dstChannel[2],Tools::_Iy());
  154.  
  155.  
  156. Mat cM = Mat::zeros(Size(image.cols, image.rows), CV_64F);
  157.  
  158. for(int i = 1; i < image.rows - 1; i++){
  159. for(int j = 1; j < image.cols - 1; j++){
  160.  
  161. g11= pow(ixR.at<double>(i,j),2)+pow(ixG.at<double>(i,j),2)+pow(ixB.at<double>(i,j),2);
  162. g12= ixR.at<double>(i,j)*iyR.at<double>(i,j)+ixG.at<double>(i,j)*iyG.at<double>(i,j)+ixB.at<double>(i,j)*iyB.at<double>(i,j);
  163. g22= pow(iyR.at<double>(i,j),2)+pow(iyG.at<double>(i,j),2)+pow(iyB.at<double>(i,j),2);
  164.  
  165. lambdaPlus=( g11 + g22 + sqrt( pow( (g11-g22),2 )+4*pow( g12,2 ) ) )/2;
  166. lambdaMin=( g11 + g22 - sqrt( pow( (g11-g22),2 )+4*pow( g12,2 ) ) )/2;
  167.  
  168. cM.at<double>(i,j) = sqrt(lambdaPlus+lambdaMin);
  169.  
  170. }
  171. }
  172.  
  173. cM.convertTo(cM, CV_8UC1, 255.0/1.0);
  174. return cM;
  175. }
  176.  
  177.  
  178. vector<Mat> Tools::_DiZenzoRestoration(Mat image,int iterations,double dt){
  179.  
  180. vector<Mat> dst(iterations+1);
  181. Mat dstChannel[3];
  182. image.copyTo(dst.at(0));
  183.  
  184.  
  185. double A,D;
  186. Mat ixyB,ixxB,iyyB,ixB,iyB;
  187. Mat ixyG,ixxG,iyyG,ixG,iyG;
  188. Mat ixyR,ixxR,iyyR,ixR,iyR;
  189. double sum,ixx,iyy,ix,iy,ixy,c;
  190.  
  191. Mat N;
  192.  
  193. for(int t=0; t<iterations;t++){cout<<t<<endl;
  194. split(dst.at(t),dstChannel);
  195.  
  196. N=_DiZenzoGradient(dst.at(t));
  197. N.convertTo(N,CV_64F,1.0/255.0);
  198.  
  199. ixyB = Tools::_Convolution(dstChannel[0],Tools::_Ixy());
  200. ixxB = Tools::_Convolution(dstChannel[0],Tools::_Ixx());
  201. iyyB = Tools::_Convolution(dstChannel[0],Tools::_Iyy());
  202. ixB = Tools::_Convolution(dstChannel[0],Tools::_Ix());
  203. iyB = Tools::_Convolution(dstChannel[0],Tools::_Iy());
  204.  
  205. ixyG = Tools::_Convolution(dstChannel[1],Tools::_Ixy());
  206. ixxG = Tools::_Convolution(dstChannel[1],Tools::_Ixx());
  207. iyyG = Tools::_Convolution(dstChannel[1],Tools::_Iyy());
  208. ixG = Tools::_Convolution(dstChannel[1],Tools::_Ix());
  209. iyG = Tools::_Convolution(dstChannel[1],Tools::_Iy());
  210.  
  211. ixyR = Tools::_Convolution(dstChannel[2],Tools::_Ixy());
  212. ixxR = Tools::_Convolution(dstChannel[2],Tools::_Ixx());
  213. iyyR = Tools::_Convolution(dstChannel[2],Tools::_Iyy());
  214. ixR = Tools::_Convolution(dstChannel[2],Tools::_Ix());
  215. iyR = Tools::_Convolution(dstChannel[2],Tools::_Iy());
  216.  
  217. dstChannel[0].convertTo(dstChannel[0],CV_64F,1.0/255.0);
  218. dstChannel[1].convertTo(dstChannel[1],CV_64F,1.0/255.0);
  219. dstChannel[2].convertTo(dstChannel[2],CV_64F,1.0/255.0);
  220.  
  221. for(int i = 1; i < image.rows-1; i++){
  222. for(int j = 1; j < image.cols-1; j++){
  223.  
  224. ixy = ixyB.at<double>(i,j);
  225. ixx = ixxB.at<double>(i,j);
  226. iyy = iyyB.at<double>(i,j);
  227. ix = ixB.at<double>(i,j);
  228. iy = iyB.at<double>(i,j);
  229.  
  230. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  231. D = (ix*ix)+(iy*iy);
  232. c = A/D;
  233.  
  234. dstChannel[0].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  235.  
  236. ixy = ixyG.at<double>(i,j);
  237. ixx = ixxG.at<double>(i,j);
  238. iyy = iyyG.at<double>(i,j);
  239. ix = ixG.at<double>(i,j);
  240. iy = iyG.at<double>(i,j);
  241.  
  242. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  243. D = (ix*ix)+(iy*iy);
  244. c = A/D;
  245.  
  246. dstChannel[1].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  247.  
  248. ixy = ixyR.at<double>(i,j);
  249. ixx = ixxR.at<double>(i,j);
  250. iyy = iyyR.at<double>(i,j);
  251. ix = ixR.at<double>(i,j);
  252. iy = iyR.at<double>(i,j);
  253.  
  254. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  255. D = (ix*ix)+(iy*iy);
  256. c = A/D;
  257.  
  258. dstChannel[2].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  259. }
  260. }
  261.  
  262. dstChannel[0].convertTo(dstChannel[0],CV_8UC1,255.0);
  263. dstChannel[1].convertTo(dstChannel[1],CV_8UC1,255.0);
  264. dstChannel[2].convertTo(dstChannel[2],CV_8UC1,255.0);
  265.  
  266. merge(dstChannel,3,dst.at(t));
  267.  
  268. dst.at(t).copyTo(dst.at(t+1));
  269. }
  270.  
  271. return dst;
  272. }
  273.  
  274. vector<Mat> Tools::_SapiroRestoration(Mat image,int iterations,double dt){
  275.  
  276. vector<Mat> dst(iterations+1);
  277. Mat dstChannel[3];
  278. image.copyTo(dst.at(0));
  279.  
  280.  
  281. double A,D;
  282. Mat ixyB,ixxB,iyyB,ixB,iyB;
  283. Mat ixyG,ixxG,iyyG,ixG,iyG;
  284. Mat ixyR,ixxR,iyyR,ixR,iyR;
  285. double sum,ixx,iyy,ix,iy,ixy,c;
  286.  
  287. Mat N;
  288.  
  289. for(int t=0; t<iterations;t++){cout<<t<<endl;
  290. split(dst.at(t),dstChannel);
  291.  
  292. N=_SapiroGradient(dst.at(t));
  293. N.convertTo(N,CV_64F,1.0/255.0);
  294.  
  295. ixyB = Tools::_Convolution(dstChannel[0],Tools::_Ixy());
  296. ixxB = Tools::_Convolution(dstChannel[0],Tools::_Ixx());
  297. iyyB = Tools::_Convolution(dstChannel[0],Tools::_Iyy());
  298. ixB = Tools::_Convolution(dstChannel[0],Tools::_Ix());
  299. iyB = Tools::_Convolution(dstChannel[0],Tools::_Iy());
  300.  
  301. ixyG = Tools::_Convolution(dstChannel[1],Tools::_Ixy());
  302. ixxG = Tools::_Convolution(dstChannel[1],Tools::_Ixx());
  303. iyyG = Tools::_Convolution(dstChannel[1],Tools::_Iyy());
  304. ixG = Tools::_Convolution(dstChannel[1],Tools::_Ix());
  305. iyG = Tools::_Convolution(dstChannel[1],Tools::_Iy());
  306.  
  307. ixyR = Tools::_Convolution(dstChannel[2],Tools::_Ixy());
  308. ixxR = Tools::_Convolution(dstChannel[2],Tools::_Ixx());
  309. iyyR = Tools::_Convolution(dstChannel[2],Tools::_Iyy());
  310. ixR = Tools::_Convolution(dstChannel[2],Tools::_Ix());
  311. iyR = Tools::_Convolution(dstChannel[2],Tools::_Iy());
  312.  
  313. dstChannel[0].convertTo(dstChannel[0],CV_64F,1.0/255.0);
  314. dstChannel[1].convertTo(dstChannel[1],CV_64F,1.0/255.0);
  315. dstChannel[2].convertTo(dstChannel[2],CV_64F,1.0/255.0);
  316.  
  317. for(int i = 1; i < image.rows-1; i++){
  318. for(int j = 1; j < image.cols-1; j++){
  319.  
  320. ixy = ixyB.at<double>(i,j);
  321. ixx = ixxB.at<double>(i,j);
  322. iyy = iyyB.at<double>(i,j);
  323. ix = ixB.at<double>(i,j);
  324. iy = iyB.at<double>(i,j);
  325.  
  326. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  327. D = (ix*ix)+(iy*iy);
  328. c = A/D;
  329.  
  330. dstChannel[0].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  331.  
  332. ixy = ixyG.at<double>(i,j);
  333. ixx = ixxG.at<double>(i,j);
  334. iyy = iyyG.at<double>(i,j);
  335. ix = ixG.at<double>(i,j);
  336. iy = iyG.at<double>(i,j);
  337.  
  338. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  339. D = (ix*ix)+(iy*iy);
  340. c = A/D;
  341.  
  342. dstChannel[1].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  343.  
  344. ixy = ixyR.at<double>(i,j);
  345. ixx = ixxR.at<double>(i,j);
  346. iyy = iyyR.at<double>(i,j);
  347. ix = ixR.at<double>(i,j);
  348. iy = iyR.at<double>(i,j);
  349.  
  350. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  351. D = (ix*ix)+(iy*iy);
  352. c = A/D;
  353.  
  354. dstChannel[2].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  355. }
  356. }
  357.  
  358. dstChannel[0].convertTo(dstChannel[0],CV_8UC1,255.0);
  359. dstChannel[1].convertTo(dstChannel[1],CV_8UC1,255.0);
  360. dstChannel[2].convertTo(dstChannel[2],CV_8UC1,255.0);
  361.  
  362. merge(dstChannel,3,dst.at(t));
  363.  
  364. dst.at(t).copyTo(dst.at(t+1));
  365. }
  366.  
  367. return dst;
  368. }
  369.  
  370. vector<Mat> Tools::_BlomgrenRestoration(Mat image,int iterations,double dt){
  371.  
  372. vector<Mat> dst(iterations+1);
  373. Mat dstChannel[3];
  374. image.copyTo(dst.at(0));
  375.  
  376.  
  377. double A,D;
  378. Mat ixyB,ixxB,iyyB,ixB,iyB;
  379. Mat ixyG,ixxG,iyyG,ixG,iyG;
  380. Mat ixyR,ixxR,iyyR,ixR,iyR;
  381. double sum,ixx,iyy,ix,iy,ixy,c;
  382.  
  383. Mat N;
  384.  
  385. for(int t=0; t<iterations;t++){cout<<t<<endl;
  386. split(dst.at(t),dstChannel);
  387.  
  388. N=_BlomgrenGradient(dst.at(t));
  389. N.convertTo(N,CV_64F,1.0/255.0);
  390.  
  391. ixyB = Tools::_Convolution(dstChannel[0],Tools::_Ixy());
  392. ixxB = Tools::_Convolution(dstChannel[0],Tools::_Ixx());
  393. iyyB = Tools::_Convolution(dstChannel[0],Tools::_Iyy());
  394. ixB = Tools::_Convolution(dstChannel[0],Tools::_Ix());
  395. iyB = Tools::_Convolution(dstChannel[0],Tools::_Iy());
  396.  
  397. ixyG = Tools::_Convolution(dstChannel[1],Tools::_Ixy());
  398. ixxG = Tools::_Convolution(dstChannel[1],Tools::_Ixx());
  399. iyyG = Tools::_Convolution(dstChannel[1],Tools::_Iyy());
  400. ixG = Tools::_Convolution(dstChannel[1],Tools::_Ix());
  401. iyG = Tools::_Convolution(dstChannel[1],Tools::_Iy());
  402.  
  403. ixyR = Tools::_Convolution(dstChannel[2],Tools::_Ixy());
  404. ixxR = Tools::_Convolution(dstChannel[2],Tools::_Ixx());
  405. iyyR = Tools::_Convolution(dstChannel[2],Tools::_Iyy());
  406. ixR = Tools::_Convolution(dstChannel[2],Tools::_Ix());
  407. iyR = Tools::_Convolution(dstChannel[2],Tools::_Iy());
  408.  
  409. dstChannel[0].convertTo(dstChannel[0],CV_64F,1.0/255.0);
  410. dstChannel[1].convertTo(dstChannel[1],CV_64F,1.0/255.0);
  411. dstChannel[2].convertTo(dstChannel[2],CV_64F,1.0/255.0);
  412.  
  413. for(int i = 1; i < image.rows-1; i++){
  414. for(int j = 1; j < image.cols-1; j++){
  415.  
  416. ixy = ixyB.at<double>(i,j);
  417. ixx = ixxB.at<double>(i,j);
  418. iyy = iyyB.at<double>(i,j);
  419. ix = ixB.at<double>(i,j);
  420. iy = iyB.at<double>(i,j);
  421.  
  422. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  423. D = (ix*ix)+(iy*iy);
  424. c = A/D;
  425.  
  426. dstChannel[0].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  427.  
  428. ixy = ixyG.at<double>(i,j);
  429. ixx = ixxG.at<double>(i,j);
  430. iyy = iyyG.at<double>(i,j);
  431. ix = ixG.at<double>(i,j);
  432. iy = iyG.at<double>(i,j);
  433.  
  434. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  435. D = (ix*ix)+(iy*iy);
  436. c = A/D;
  437.  
  438. dstChannel[1].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  439.  
  440. ixy = ixyR.at<double>(i,j);
  441. ixx = ixxR.at<double>(i,j);
  442. iyy = iyyR.at<double>(i,j);
  443. ix = ixR.at<double>(i,j);
  444. iy = iyR.at<double>(i,j);
  445.  
  446. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  447. D = (ix*ix)+(iy*iy);
  448. c = A/D;
  449.  
  450. dstChannel[2].at<double>(i,j) += exp(-N.at<double>(i,j))*(dt*c);
  451. }
  452. }
  453.  
  454. dstChannel[0].convertTo(dstChannel[0],CV_8UC1,255.0);
  455. dstChannel[1].convertTo(dstChannel[1],CV_8UC1,255.0);
  456. dstChannel[2].convertTo(dstChannel[2],CV_8UC1,255.0);
  457.  
  458. merge(dstChannel,3,dst.at(t));
  459.  
  460. dst.at(t).copyTo(dst.at(t+1));
  461. }
  462.  
  463. return dst;
  464. }
  465.  
  466. vector<Mat> Tools::_DericheRestoration(Mat image,int iterations,double dt){
  467.  
  468. vector<Mat> dst(iterations+1);
  469. Mat dstChannel[3];
  470. image.copyTo(dst.at(0));
  471.  
  472. double A,D;
  473. Mat ixyB,ixxB,iyyB,ixB,iyB;
  474. Mat ixyG,ixxG,iyyG,ixG,iyG;
  475. Mat ixyR,ixxR,iyyR,ixR,iyR;
  476.  
  477. double sum,ixx,iyy,ix,iy,ixy,ee,nn;
  478.  
  479. Mat N;
  480.  
  481. for(int t=0; t<iterations;t++){cout<<t<<endl;
  482. split(dst.at(t),dstChannel);
  483.  
  484. N=_DiZenzoGradient(dst.at(t));
  485. N.convertTo(N,CV_64F,1.0/255.0);
  486.  
  487. ixyB = Tools::_Convolution(dstChannel[0],Tools::_Ixy());
  488. ixxB = Tools::_Convolution(dstChannel[0],Tools::_Ixx());
  489. iyyB = Tools::_Convolution(dstChannel[0],Tools::_Iyy());
  490. ixB = Tools::_Convolution(dstChannel[0],Tools::_Ix());
  491. iyB = Tools::_Convolution(dstChannel[0],Tools::_Iy());
  492.  
  493. ixyG = Tools::_Convolution(dstChannel[1],Tools::_Ixy());
  494. ixxG = Tools::_Convolution(dstChannel[1],Tools::_Ixx());
  495. iyyG = Tools::_Convolution(dstChannel[1],Tools::_Iyy());
  496. ixG = Tools::_Convolution(dstChannel[1],Tools::_Ix());
  497. iyG = Tools::_Convolution(dstChannel[1],Tools::_Iy());
  498.  
  499. ixyR = Tools::_Convolution(dstChannel[2],Tools::_Ixy());
  500. ixxR = Tools::_Convolution(dstChannel[2],Tools::_Ixx());
  501. iyyR = Tools::_Convolution(dstChannel[2],Tools::_Iyy());
  502. ixR = Tools::_Convolution(dstChannel[2],Tools::_Ix());
  503. iyR = Tools::_Convolution(dstChannel[2],Tools::_Iy());
  504.  
  505. dstChannel[0].convertTo(dstChannel[0],CV_64F,1.0/255.0);
  506. dstChannel[1].convertTo(dstChannel[1],CV_64F,1.0/255.0);
  507. dstChannel[2].convertTo(dstChannel[2],CV_64F,1.0/255.0);
  508.  
  509. for(int i = 1; i < image.rows-1; i++){
  510. for(int j = 1; j < image.cols-1; j++){
  511.  
  512. ixy = ixyB.at<double>(i,j);
  513. ixx = ixxB.at<double>(i,j);
  514. iyy = iyyB.at<double>(i,j);
  515. ix = ixB.at<double>(i,j);
  516. iy = iyB.at<double>(i,j);
  517.  
  518. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  519. D = (ix*ix)+(iy*iy);
  520. ee = A/D;
  521. A = (ixx*(ix*ix)) + (2*ixy*ix*iy) + (iyy*(iy*iy));
  522. nn = A/D;
  523.  
  524. dstChannel[0].at<double>(i,j) += (exp(-N.at<double>(i,j))*(dt*nn)+(dt*ee));
  525.  
  526. ixy = ixyG.at<double>(i,j);
  527. ixx = ixxG.at<double>(i,j);
  528. iyy = iyyG.at<double>(i,j);
  529. ix = ixG.at<double>(i,j);
  530. iy = iyG.at<double>(i,j);
  531.  
  532. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  533. D = (ix*ix)+(iy*iy);
  534. ee = A/D;
  535. A = (ixx*(ix*ix)) + (2*ixy*ix*iy) + (iyy*(iy*iy));
  536. nn = A/D;
  537.  
  538. dstChannel[1].at<double>(i,j) += (exp(-N.at<double>(i,j))*(dt*nn)+(dt*ee));
  539.  
  540. ixy = ixyR.at<double>(i,j);
  541. ixx = ixxR.at<double>(i,j);
  542. iyy = iyyR.at<double>(i,j);
  543. ix = ixR.at<double>(i,j);
  544. iy = iyR.at<double>(i,j);
  545.  
  546. A = (ixx*(iy*iy)) - (2*ixy*ix*iy) + (iyy*(ix*ix));
  547. D = (ix*ix)+(iy*iy);
  548. ee = A/D;
  549. A = (ixx*(ix*ix)) + (2*ixy*ix*iy) + (iyy*(iy*iy));
  550. nn = A/D;
  551.  
  552. dstChannel[2].at<double>(i,j) += (exp(-N.at<double>(i,j))*(dt*nn)+(dt*ee));
  553. }
  554. }
  555.  
  556. dstChannel[0].convertTo(dstChannel[0],CV_8UC1,255.0);
  557. dstChannel[1].convertTo(dstChannel[1],CV_8UC1,255.0);
  558. dstChannel[2].convertTo(dstChannel[2],CV_8UC1,255.0);
  559.  
  560. merge(dstChannel,3,dst.at(t));
  561.  
  562. dst.at(t).copyTo(dst.at(t+1));
  563. }
  564.  
  565. return dst;
  566. }
  567.  
  568.  
  569.  
  570. vector<Mat> Tools::_BelfkihRestoration(Mat image,int iterations,double dt){
  571.  
  572.  
  573. vector<Mat> dst(iterations+1);
  574. Mat dstChannel[3];
  575. image.copyTo(dst.at(0));
  576.  
  577. double A,D;
  578. Mat _ixyB,_ixxB,_iyyB,_ixB,_iyB;
  579. Mat _ixyG,_ixxG,_iyyG,_ixG,_iyG;
  580. Mat _ixyR,_ixxR,_iyyR,_ixR,_iyR;
  581.  
  582. double sum,ixx,iyy,ix,iy,ixy,ee,nn;
  583.  
  584. double ixyB,ixxB,iyyB,ixB,iyB;
  585. double ixyG,ixxG,iyyG,ixG,iyG;
  586. double ixyR,ixxR,iyyR,ixR,iyR;
  587.  
  588.  
  589. double g11,g12,g21,g22;
  590. double lambdaMin,lambdaMax;
  591. double xG2,xB2,xR2,yG2,yB2,yR2;
  592. double dteta_x,dteta_y,dx_s,dy_s,A2,lap;
  593. double dx,dy,Kval=0.0;
  594. double eeB,eeG,eeR;
  595. double nnB,nnG,nnR;
  596.  
  597. Mat _GB,_GG,_GR;
  598. Mat _LBX,_LGX,_LRX;
  599. Mat _LBY,_LGY,_LRY;
  600.  
  601. Mat N;
  602. Mat cMB = Mat::zeros(Size(image.cols, image.rows), CV_64F);
  603. Mat cMG = Mat::zeros(Size(image.cols, image.rows), CV_64F);
  604. Mat cMR = Mat::zeros(Size(image.cols, image.rows), CV_64F);
  605.  
  606. for(int t=0; t<iterations;t++){cout<<t<<endl;
  607. split(dst.at(t),dstChannel);
  608.  
  609. N=_DiZenzoGradient(dst.at(t));
  610. N.convertTo(N,CV_64F,1.0/255.0);
  611.  
  612. _ixyB = Tools::_Convolution(dstChannel[0],Tools::_Ixy());
  613. _ixxB = Tools::_Convolution(dstChannel[0],Tools::_Ixx());
  614. _iyyB = Tools::_Convolution(dstChannel[0],Tools::_Iyy());
  615. _ixB = Tools::_Convolution(dstChannel[0],Tools::_Ix());
  616. _iyB = Tools::_Convolution(dstChannel[0],Tools::_Iy());
  617.  
  618. _ixyG = Tools::_Convolution(dstChannel[1],Tools::_Ixy());
  619. _ixxG = Tools::_Convolution(dstChannel[1],Tools::_Ixx());
  620. _iyyG = Tools::_Convolution(dstChannel[1],Tools::_Iyy());
  621. _ixG = Tools::_Convolution(dstChannel[1],Tools::_Ix());
  622. _iyG = Tools::_Convolution(dstChannel[1],Tools::_Iy());
  623.  
  624. _ixyR = Tools::_Convolution(dstChannel[2],Tools::_Ixy());
  625. _ixxR = Tools::_Convolution(dstChannel[2],Tools::_Ixx());
  626. _iyyR = Tools::_Convolution(dstChannel[2],Tools::_Iyy());
  627. _ixR = Tools::_Convolution(dstChannel[2],Tools::_Ix());
  628. _iyR = Tools::_Convolution(dstChannel[2],Tools::_Iy());
  629.  
  630. _GB=Tools::_Convolution(dstChannel[0],Tools::_GaussianKernel(3,1.5));
  631. _GG=Tools::_Convolution(dstChannel[1],Tools::_GaussianKernel(3,1.5));
  632. _GR=Tools::_Convolution(dstChannel[2],Tools::_GaussianKernel(3,1.5));
  633.  
  634. _LBX=Tools::_Convolution(_GB,Tools::_Ixx());
  635. _LBY=Tools::_Convolution(_GB,Tools::_Iyy());
  636.  
  637. _LGX=Tools::_Convolution(_GG,Tools::_Ixx());
  638. _LGY=Tools::_Convolution(_GG,Tools::_Iyy());
  639.  
  640. _LRX=Tools::_Convolution(_GR,Tools::_Ixx());
  641. _LRY=Tools::_Convolution(_GR,Tools::_Iyy());
  642.  
  643. for(int i=0;i<image.rows;i++){
  644. for(int j=0;j<image.cols;j++){
  645. cMB.at<double>(i,j)=sqrt(pow(_LBX.at<double>(i,j),2)+pow(_LBX.at<double>(i,j),2));
  646. cMG.at<double>(i,j)=sqrt(pow(_LGX.at<double>(i,j),2)+pow(_LGX.at<double>(i,j),2));
  647. cMR.at<double>(i,j)=sqrt(pow(_LRX.at<double>(i,j),2)+pow(_LRX.at<double>(i,j),2));
  648. }
  649. }
  650. dstChannel[0].convertTo(dstChannel[0],CV_64F,1.0/255.0);
  651. dstChannel[1].convertTo(dstChannel[1],CV_64F,1.0/255.0);
  652. dstChannel[2].convertTo(dstChannel[2],CV_64F,1.0/255.0);
  653.  
  654. for(int i = 1; i < image.rows-1; i++){
  655. for(int j = 1; j < image.cols-1; j++){
  656.  
  657. ixyB = _ixyB.at<double>(i,j);
  658. ixxB = _ixxB.at<double>(i,j);
  659. iyyB = _iyyB.at<double>(i,j);
  660. ixB = _ixB.at<double>(i,j);
  661. iyB = _iyB.at<double>(i,j);
  662.  
  663. A = (ixxB*(iyB*iyB)) - (2*ixyB*ixB*iyB) + (iyyB*(ixB*ixB));
  664. D = (ixB*ixB)+(iyB*iyB);
  665. eeB = A/D;
  666. A = (ixxB*(ixB*ixB)) + (2*ixyB*ixB*iyB) + (iyyB*(iyB*iyB));
  667. nnB = A/D;
  668.  
  669. //dstChannel[0].at<double>(i,j) += (exp(-N.at<double>(i,j))*(dt*nnB)+(dt*eeB));
  670.  
  671. ixyG = _ixyG.at<double>(i,j);
  672. ixxG = _ixxG.at<double>(i,j);
  673. iyyG = _iyyG.at<double>(i,j);
  674. ixG = _ixG.at<double>(i,j);
  675. iyG = _iyG.at<double>(i,j);
  676.  
  677. A = (ixxG*(iyG*iyG)) - (2*ixyG*ixG*iyG) + (iyyG*(ixG*ixG));
  678. D = (ixG*ixG)+(iyG*iyG);
  679. eeG = A/D;
  680. A = (ixxG*(ixG*ixG)) + (2*ixyG*ixG*iyG) + (iyyG*(iyG*iyG));
  681. nnG = A/D;
  682.  
  683. //dstChannel[1].at<double>(i,j) += (exp(-N.at<double>(i,j))*(dt*nnG)+(dt*eeG));
  684.  
  685. ixyR = _ixyR.at<double>(i,j);
  686. ixxR = _ixxR.at<double>(i,j);
  687. iyyR = _iyyR.at<double>(i,j);
  688. ixR = _ixR.at<double>(i,j);
  689. iyR = _iyR.at<double>(i,j);
  690.  
  691. A = (ixxR*(iyR*iyR)) - (2*ixyR*ixR*iyR) + (iyyR*(ixR*ixR));
  692. D = (ixR*ixR)+(iyR*iyR);
  693. eeR = A/D;
  694. A = (ixxR*(ixR*ixR)) + (2*ixyR*ixR*iyR) + (iyyR*(iyR*iyR));
  695. nnR = A/D;
  696.  
  697. //dstChannel[2].at<double>(i,j) += (exp(-N.at<double>(i,j))*(dt*nnR)+(dt*eeR));
  698. g11=ixR*ixR+ixG*ixG+ixB*ixB;
  699. g22=iyR*iyR+iyG*iyG+iyB*iyB;
  700. g12=ixR*iyR+ixG*iyG+ixB*iyB;
  701.  
  702. A2=A*A;
  703. lap=A2 + 4*g12*g12;
  704. if(lap==0)lap=1;
  705. lambdaMax=(g11 + g22 + sqrt(A2 + 4*g12*g12))/2;
  706.  
  707. A=ixB*ixB - iyB*iyB + ixR*ixR - iyR*iyR + ixG*ixG - iyG*iyG;
  708.  
  709. dteta_x=(A*ixB*ixyB + A*ixxB*iyB - 2*ixG*ixyG*pow(iyG,2) + 2*pow(ixB,2)*ixxB*iyB +
  710. 2*ixB*ixyB*pow(iyB,2) - 2*ixB*iyB*ixR*ixxR + A*ixR*ixyR - 2*ixB*ixxB*ixR*iyR +
  711. 2*ixyB*iyB*ixR*iyR + A*ixxR*iyR - 2*pow(ixR,2)*ixxR*iyR + 2*ixB*iyB*ixyR*iyR +
  712. 2*ixR*ixyR*pow(iyR,2) - 2*ixB*iyB*ixG*ixxG - 2*ixR*iyR*ixG*ixxG + A*ixG*ixyG - 2*ixB*ixxB*ixG*iyG +
  713. 2*ixR*ixxR*ixG*iyG + 2*ixyB*iyB*ixG*iyG - A*Xxx*iyG - 2*pow(ixG,2)*ixxG*iyG + 2*ixyR*iyR*ixG*iyG +
  714. 2*ixB*iyB*ixyG*iyG + 2*ixR*iyR*ixyG*iyG) / (pow(A,2) + 4*pow(ixB*iyB + ixR*iyR + ixG*iyG, 2));
  715.  
  716. dteta_y=(A*ixyB*iyB - 2*ixB*ixyB*iyB + A*ixB*iyyB + 2*ixB*iyB*iyB*iyyB - 2*ixB*iyB*ixR*ixyR - 2*ixB*ixyB*ixR*iyR + //??
  717. 2*iyB*iyyB*ixR*iyR + A*ixyR*iyB - 2*ixR*ixR*ixyR*iyR + 2*ixB*iyB*iyR*iyyR + 2*ixR*iyR*iyR*iyyR - 2*ixB*iyB*ixG*iyG -
  718. 2*ixR*iyR*ixG*ixyG - 2*ixB*ixyB*ixG*iyG + 2*iyB*iyyB*ixG*iyG - 2*ixR*ixyR*ixyR*ixG*iyG + 2*iyR*iyyR*ixG*iyG + A*ixyG*iyG -
  719. 2*ixG*ixG*ixyG*iyG + A*ixG*iyyG + 2*ixB*iyB*iyG*iyyG + 2*ixR*iyR*ixG*iyyG + 2*ixG*iyG*iyG*iyyG + //??
  720. 2*ixB*ixyB*iyB - 2*ixB*ixyB*iyB) / (pow(A,2) + 4*pow(ixB*iyB + ixR*iyR + ixG*iyG, 2));
  721.  
  722. dx_s=1/(sqrt(pow(g12,2) + pow(g11-lambdaMax,2)) * (g11 - lambdaMax));
  723.  
  724. dy_s=1/(sqrt(pow(g12,2) + pow(g11*lambdaMax,2)) * g12);
  725.  
  726. Kval = dteta_x * dx_s + dteta_y * dy_s;
  727.  
  728. dstChannel[0].at<double>(i,j) += (gF(cMB.at<double>(i,j)) * (nnB) +gF(Kval)* (eeB))*dt;
  729. dstChannel[1].at<double>(i,j) += (gF(cMG.at<double>(i,j)) * (nnG) +gF(Kval)* (eeG))*dt;
  730. dstChannel[2].at<double>(i,j) += (gF(cMR.at<double>(i,j)) * (nnR) +gF(Kval)* (eeR))*dt;
  731. }
  732. }
  733.  
  734. dstChannel[0].convertTo(dstChannel[0],CV_8UC1,255.0);
  735. dstChannel[1].convertTo(dstChannel[1],CV_8UC1,255.0);
  736. dstChannel[2].convertTo(dstChannel[2],CV_8UC1,255.0);
  737.  
  738. merge(dstChannel,3,dst.at(t));
  739.  
  740. dst.at(t).copyTo(dst.at(t+1));
  741. }
  742.  
  743. return dst;
  744. }
  745.  
  746.  
  747. double Tools::gaussianOp(int x,int y,double sigma){
  748. return ( (1/(2*M_PI*pow(sigma,2))) * exp(-(pow(x,2)+pow(y,2)/(2*pow(sigma,2))))) ;
  749. }
  750.  
  751. Mat Tools::_GaussianKernel(int size,double sigma){
  752. Mat gkernel= Mat(size,size,CV_64F);
  753. double sumNorm=0;
  754. for(int i=0;i<size;i++){
  755. for(int j=0;j<size;j++){
  756. gkernel.at<double>(i,j)=gaussianOp(i-size/2,j-size/2,sigma);
  757. sumNorm+=gkernel.at<double>(i,j);
  758. }
  759. }
  760.  
  761. for(int i=0;i<size;i++){
  762. for(int j=0;j<size;j++){
  763. gkernel.at<double>(i,j)/=sumNorm;
  764. }
  765. }
  766.  
  767. return gkernel;
  768. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!