Guest User

Untitled

a guest
Sep 25th, 2020
32
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. #include "mainwindow.h"
  2.  
  3. using namespace std;
  4. using namespace cv;
  5.  
  6. // Найти уравнение прямой
  7. void line_equation(Point src1, Point src2, double& k, int& b) {
  8. if ((src1.x != src2.x) && (src1.y != src2.y)) {
  9. k = (double) (src1.y - src2.y) / (double) (src1.x - src2.x);
  10. b = src1.y - (k * src1.x);
  11. } else {
  12. k = 0;
  13. }
  14. }
  15.  
  16. void find_parallel(Point src, vector<Point>& dst, double k, int rows) {
  17. int b = src.y - k * src.x;
  18. dst.push_back(Point(-b/k, 0));
  19. dst.push_back(Point((rows - b) / k, rows));
  20. }
  21.  
  22. // Найти точки пересечения прямой и контура
  23. void line_interseption(vector<Point> src_line, Point p, Point q, vector<Point>& dst_point) {
  24. int d = 0;
  25. bool first_point_flag = false;
  26. bool second_point_flag = false;
  27. for(size_t i = 1; i < src_line.size(); i++) {
  28. d = abs((p.y - q.y)*src_line[i].x - (p.x - q.x)*src_line[i].y + p.x * q.y - p.y * q.x)
  29. /sqrt((p.y - q.y) * (p.y - q.y) + (p.x - q.x) * (p.x - q.x));
  30. if(p.x != q.x) {
  31. if((d == 0) && !first_point_flag && (p.x > src_line[i].x)) {
  32. dst_point.push_back(src_line[i]);
  33. first_point_flag = true;
  34. }
  35. if((d == 0) && !second_point_flag && (p.x < src_line[i].x)) {
  36. dst_point.push_back(src_line[i]);
  37. second_point_flag = true;
  38. }
  39. } else {
  40. if((d == 0) && !first_point_flag && (p.y > src_line[i].y)) {
  41. dst_point.push_back(src_line[i]);
  42. first_point_flag = true;
  43. }
  44. if((d == 0) && !second_point_flag && (p.y < src_line[i].y)) {
  45. dst_point.push_back(src_line[i]);
  46. second_point_flag = true;
  47. }
  48. }
  49. }
  50. }
  51.  
  52. // Найти точку пересечения прямых
  53. void intersection_points(Point flp1, Point flp2, Point slp1, Point slp2, Point& dst_point) {
  54. double k1 = (double)(flp2.y - flp1.y) / (flp2.x - flp1.x);
  55. double k2 = (double)(slp2.y - slp1.y) / (slp2.x - slp1.x);
  56. double b1 = k1 * flp1.x - flp1.y;
  57. double b2 = k2 * slp1.x - slp1.y;
  58. dst_point.x = -(b1 - b2) / (k2 - k1);
  59. dst_point.y = -(k2 * b1 - k1 * b2) / (k2 - k1);
  60. }
  61.  
  62.  
  63. // Найти максимально удаленную точку от прямой
  64. void find_max(vector<Point> sepCont, Point p, Point q, Point &dst) {
  65. dst = sepCont[0];
  66. int d = 0;
  67. int d_max = abs((p.y - q.y)*sepCont[0].x - (p.x - q.x)*sepCont[0].y + p.x * q.y - p.y * q.x)
  68. /sqrt((p.y - q.y) * (p.y - q.y) + (p.x - q.x) * (p.x - q.x));
  69. for(size_t i = 1; i < sepCont.size(); i++) {
  70. d = abs((p.y - q.y)*sepCont[i].x - (p.x - q.x)*sepCont[i].y + p.x * q.y - p.y * q.x)
  71. /sqrt((p.y - q.y) * (p.y - q.y) + (p.x - q.x) * (p.x - q.x));
  72. if(d > d_max) {
  73. dst = sepCont[i];
  74. d_max = d;
  75. }
  76. }
  77. }
  78.  
  79. // Найти расстаяние между двумя точками
  80. int find_dist(Point p1, Point p2) {
  81. return sqrt((p2.x - p1.x) * (p2.x - p1.x) + (p2.y - p1.y) * (p2.y - p1.y));
  82. }
  83.  
  84. // Найти точку на прямой
  85. void point_on_line(Point src1, Point src2, double k, int b, Point& dst) {
  86. double dif = 0;
  87. if(src1.y - src2.y != 0) {
  88. dif = (double)abs(src1.x - src2.x) / (double)abs(src1.y - src2.y);
  89. } else {
  90. dif = 2;
  91. }
  92. int dist = find_dist(src1, src2);
  93. if(dif > 1) {
  94. if(src1.x > src2.x) {
  95. dst.y = k * (src1.x - dist/10) + b;
  96. dst.x = src1.x - dist/10;
  97. } else {
  98. dst.y = k * (src1.x + dist/10) + b;
  99. dst.x = src1.x + dist/10;
  100. }
  101. } else {
  102. if(src1.y > src2.y) {
  103. dst.y = src1.y - dist/10;
  104. dst.x = (dst.y - b) / k;
  105. } else {
  106. dst.y = src1.y + dist/10;
  107. dst.x = (dst.y - b) / k;
  108. }
  109. }
  110. }
  111.  
  112. void point_on_line(vector<Point> src, double k, int b, vector<Point>& dst) {
  113. double dif = 0;
  114. if(src[0].y - src[1].y != 0) {
  115. dif = (double)abs(src[0].x - src[1].x) / (double)abs(src[0].y - src[1].y);
  116. } else {
  117. dif = 2;
  118. }
  119. int dist = find_dist(src[0], src[1]);
  120. if(dif > 1) {
  121. if(src[0].x > src[1].x) {
  122. dst[0].y = k * (src[0].x - dist/10) + b;
  123. dst[0].x = src[0].x - dist/10;
  124. dst[1].y = k * (src[1].x + dist/10) + b;
  125. dst[1].x = src[1].x + dist/10;
  126. } else {
  127. dst[0].y = k * (src[0].x + dist/10) + b;
  128. dst[0].x = src[0].x + dist/10;
  129. dst[1].y = k * (src[1].x - dist/10) + b;
  130. dst[1].x = src[1].x - dist/10;
  131. }
  132. } else {
  133. if(src[0].y > src[1].y) {
  134. dst[0].y = src[0].y - dist/10;
  135. dst[0].x = (dst[0].y - b) / k;
  136. dst[1].y = src[1].y + dist/10;
  137. dst[1].x = (dst[1].y - b) / k;
  138. } else {
  139. dst[0].y = src[0].y + dist/10;
  140. dst[0].x = (dst[0].y - b) / k;
  141. dst[1].y = src[1].y - dist/10;
  142. dst[1].x = (dst[1].y - b) / k;
  143. }
  144. }
  145. }
  146.  
  147. // Извлечь прямоугольник по 4 точкам
  148. void extractRect(Mat src, Mat dst, int min_y, int min_x, int max_y, int max_x,
  149. int b1, int b2, int b3, int b4, double k1, double k2) {
  150. while(min_y <= max_y) {
  151. int x1 = (min_y - b3) / k1;
  152. int x2 = (min_y - b1) / k2;
  153. if(x1 < min_x)
  154. x1 = (min_y - b2) / k1;
  155. if(x2 > max_x)
  156. x2 = (min_y - b4) / k2;
  157. if(x1 > x2) {
  158. int temp = x2;
  159. x2 = x1;
  160. x1 = temp;
  161. }
  162. for(int i = x1; i <= x2; i++) {
  163. dst.at<uchar>(Point(i, min_y)) = src.at<uchar>(Point(i, min_y));
  164. }
  165. min_y++;
  166. }
  167. }
  168.  
  169. void separeteContour(vector<vector<Point>> src, vector<Point>& top, vector<Point>& low, double k, int b) {
  170. for(size_t i = 0; i < src.size(); i++) {
  171. for(size_t j = 0; j < src[i].size(); j++) {
  172. int y_cont = k * src[i][j].x + b;
  173. if(src[i][j].y < y_cont) {
  174. low.push_back(src[i][j]);
  175. } else {
  176. top.push_back(src[i][j]);
  177. }
  178. }
  179. }
  180. }
  181.  
  182. int main()
  183. {
  184. Mat fg, bg, blured, preform, bin, ROI, erosion_dst, cont;
  185. vector<vector<Point>> contours, true_cont;
  186. Point2d vertices[4];
  187. int b_main, b1, b2, b3, b4;
  188. double k_main, k_peak;
  189.  
  190. bg = imread("C:\\2\\bg.bmp", IMREAD_GRAYSCALE);
  191. fg = imread("C:\\2\\pf1.bmp", IMREAD_GRAYSCALE);
  192.  
  193. resize(fg, fg, Size(), 0.3, 0.3, INTER_AREA);
  194. resize(bg, bg, Size(), 0.3, 0.3, INTER_AREA);
  195.  
  196. absdiff(fg, bg, preform);
  197.  
  198. threshold(preform, bin, 10, 255, THRESH_BINARY);
  199. findContours(bin, contours, RETR_EXTERNAL, CHAIN_APPROX_NONE);
  200.  
  201. uint32_t resolution = fg.rows * fg.cols;
  202. uint32_t kernelSize = resolution < 1280 * 1280 ? 7:
  203. resolution < 2000 * 2000 ? 15:
  204. resolution < 3000 * 3000 ? 27:
  205. 45;
  206.  
  207. GaussianBlur(bin, blured, Size(kernelSize, kernelSize), 9, 9);
  208.  
  209. imshow("bin", bin);
  210. imshow("pf", preform);
  211.  
  212. for(size_t i = 0; i < contours.size(); i++) {
  213. if(contours[i].size() > 1000)
  214. true_cont.push_back(contours[i]);
  215. }
  216.  
  217. vector<Vec3f> circles;
  218. HoughCircles(blured, circles, HOUGH_GRADIENT, 1, fg.rows/16, 100, 30, 1, 100);
  219.  
  220. Moments mnt = moments(true_cont[0]);
  221. Point center_mass(mnt.m10/mnt.m00, mnt.m01/mnt.m00);
  222. Point circle_center(cvRound(circles[0][0]), cvRound(circles[0][1]));
  223.  
  224. line_equation(center_mass, circle_center, k_main, b_main);
  225.  
  226. int main_dist = find_dist(center_mass, circle_center);
  227.  
  228. Point peak1, peak2;
  229. vector<Point> top_cont, low_cont;
  230. separeteContour(true_cont, top_cont, low_cont, k_main, b_main);
  231.  
  232. find_max(top_cont, center_mass, circle_center, peak1);
  233. find_max(low_cont, center_mass, circle_center, peak2);
  234.  
  235. line_equation(peak1, peak2, k_peak, b1);
  236.  
  237. Point neck_point;
  238.  
  239. intersection_points(center_mass, circle_center, peak1, peak2, neck_point);
  240.  
  241. point_on_line(neck_point, center_mass, k_main, b_main, neck_point);
  242.  
  243. line(fg, Point(-b_main/k_main, 0), Point((fg.rows-b_main)/k_main, fg.rows), 255, 1);
  244.  
  245. vector<Point> neck_parallel;
  246. vector<Point> circle_parallel;
  247. vector<Point> circle_int;
  248. vector<Point> circle_int_par1;
  249. vector<Point> circle_int_par2;
  250.  
  251. find_parallel(neck_point, neck_parallel, k_peak, fg.rows);
  252. find_parallel(circle_center, circle_parallel, k_peak, fg.rows);
  253.  
  254. line_interseption(top_cont, circle_parallel[0], circle_parallel[1], circle_int);
  255. line_interseption(low_cont, circle_parallel[0], circle_parallel[1], circle_int);
  256.  
  257. b2 = circle_int[0].y - k_peak * circle_int[0].x;
  258.  
  259. point_on_line(circle_int, k_peak, b2, circle_int);
  260.  
  261. find_parallel(circle_int[0], circle_int_par1, k_main, fg.rows);
  262. find_parallel(circle_int[1], circle_int_par2, k_main, fg.rows);
  263.  
  264. Point neck_int1, neck_int2;
  265. intersection_points(circle_int_par1[0], circle_int_par1[1], neck_parallel[0], neck_parallel[1], neck_int1);
  266. intersection_points(circle_int_par2[0], circle_int_par2[1], neck_parallel[0], neck_parallel[1], neck_int2);
  267.  
  268. vertices[0] = circle_int[0];
  269. vertices[1] = circle_int[1];
  270. vertices[2] = neck_int1;
  271. vertices[3] = neck_int2;
  272.  
  273. circle(fg, neck_point, 2, 0, 2);
  274. circle(fg, peak1, 2, 100, 2);
  275. circle(fg, peak2, 2, 100, 2);
  276. circle(fg, vertices[0], 2, 0, 2);
  277. circle(fg, vertices[1], 2, 0, 2);
  278. circle(fg, vertices[2], 2, 0, 2);
  279. circle(fg, vertices[3], 2, 0, 2);
  280. line(fg, neck_parallel[1], neck_parallel[0], 255, 1);
  281. line(fg, circle_parallel[1], circle_parallel[0], 255, 1);
  282. line(fg, circle_int_par1[1], circle_int_par1[0], 255, 1);
  283. line(fg, circle_int_par2[1], circle_int_par2[0], 255, 1);
  284.  
  285. imshow("fg", fg);
  286. imshow("pf", preform);
  287. imshow("bl", blured);
  288.  
  289. waitKey();
  290.  
  291. return EXIT_SUCCESS;
  292. }
  293.  
RAW Paste Data