OTSImageProcess.cpp 42 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289129012911292129312941295129612971298129913001301130213031304130513061307130813091310131113121313131413151316131713181319132013211322132313241325132613271328132913301331133213331334133513361337133813391340134113421343134413451346134713481349135013511352135313541355135613571358135913601361136213631364136513661367136813691370137113721373137413751376137713781379138013811382138313841385138613871388138913901391139213931394139513961397139813991400140114021403140414051406140714081409141014111412141314141415141614171418141914201421142214231424142514261427142814291430143114321433143414351436143714381439144014411442144314441445144614471448144914501451145214531454145514561457145814591460146114621463146414651466146714681469147014711472147314741475147614771478147914801481148214831484148514861487148814891490149114921493149414951496149714981499150015011502150315041505150615071508150915101511151215131514151515161517151815191520152115221523152415251526152715281529153015311532153315341535153615371538153915401541154215431544154515461547154815491550155115521553155415551556155715581559156015611562156315641565156615671568156915701571157215731574157515761577157815791580158115821583158415851586158715881589159015911592159315941595159615971598159916001601160216031604160516061607160816091610161116121613161416151616161716181619
  1. #pragma once
  2. #include "stdafx.h"
  3. #include <opencv2/core/core.hpp>
  4. #include <opencv2/highgui/highgui.hpp>
  5. #include <opencv2/opencv.hpp>
  6. #include "OTSImageProcess.h"
  7. #include "OTSImageProcessParam.h"
  8. #include <OTSFieldData.h>
  9. #include "OTSMorphology.h"
  10. #include "../OTSLog/COTSUtilityDllFunExport.h"
  11. #include "FieldMgr.h"
  12. #include "BaseFunction.h"
  13. namespace OTSIMGPROC
  14. {
  15. using namespace cv;
  16. using namespace std;
  17. const int nBlackColor = 255;
  18. COTSImageProcess::COTSImageProcess(COTSImageProcessParamPtr a_pImageProcessParam)
  19. {
  20. m_imageProcessParam = a_pImageProcessParam;
  21. }
  22. COTSImageProcess::~COTSImageProcess()
  23. {
  24. }
  25. BOOL COTSImageProcess::RemoveBGByCVconnectivities(CBSEImgPtr inBSEImg, double a_pixelSize, COTSFieldDataPtr m_pFieldData)
  26. {
  27. ASSERT(m_pFieldData);
  28. ASSERT(inBSEImg);
  29. ASSERT(m_imageProcessParam);
  30. int nWidthImg = inBSEImg->GetWidth();
  31. int nHeightImg = inBSEImg->GetHeight();
  32. m_pFieldData->Width = nWidthImg;
  33. m_pFieldData->Height = nHeightImg;
  34. long nImgSize = nWidthImg * nHeightImg;
  35. BYTE* pSrcImg = inBSEImg->GetImageDataPointer();
  36. BYTE* pTempImg = new BYTE[nImgSize];
  37. CRect r = CRect(0, 0, nWidthImg, nHeightImg);
  38. CBSEImgPtr imgNoBGBinary = CBSEImgPtr(new CBSEImg(r));
  39. long nNumParticle = 0;
  40. RemoveBackGround(inBSEImg, m_imageProcessParam, imgNoBGBinary, nNumParticle);
  41. BYTE* pPixel = imgNoBGBinary->GetImageDataPointer();
  42. long nPtStart = m_imageProcessParam->GetParticleGray().GetStart();
  43. long nPtEnd = m_imageProcessParam->GetParticleGray().GetEnd();
  44. if (nNumParticle == 0)
  45. {
  46. COTSParticleList listParticleEmpty;
  47. listParticleEmpty.clear();
  48. m_pFieldData->SetParticleList(listParticleEmpty);
  49. }
  50. else
  51. {
  52. // get the area image
  53. Mat cvcopyImg = Mat(nHeightImg, nWidthImg, CV_8UC1, pPixel);
  54. Mat labels = Mat::zeros(cvcopyImg.size(), CV_32S);
  55. Mat stats, centroids;
  56. int number = connectedComponentsWithStats(cvcopyImg, labels, stats, centroids, 8, CV_32S);
  57. double rMin = m_imageProcessParam->GetIncArea().GetStart()/2.0;
  58. double rMax = m_imageProcessParam->GetIncArea().GetEnd()/2.0;
  59. double partAreaMin = rMin * rMin * 3.14159;
  60. double partAreaMax = rMax * rMax * 3.14159;
  61. COTSParticleList listParticleOut;
  62. for (size_t i = 1; i < number; i++)
  63. {
  64. int center_x = centroids.at<double>(i, 0);
  65. int center_y = centroids.at<double>(i, 1);
  66. //앤近긋움
  67. int x = stats.at<int>(i, CC_STAT_LEFT);
  68. int y = stats.at<int>(i, CC_STAT_TOP);
  69. int w = stats.at<int>(i, CC_STAT_WIDTH);
  70. int h = stats.at<int>(i, CC_STAT_HEIGHT);
  71. int area = stats.at<int>(i, CC_STAT_AREA);
  72. double actualArea = area * a_pixelSize * a_pixelSize;
  73. if (actualArea >= partAreaMin && actualArea < partAreaMax)
  74. {
  75. Rect rectMax = Rect(x, y, w, h);
  76. Mat rectROI = labels(rectMax).clone();
  77. Mat imageROI = Mat::zeros(rectMax.size(), cvcopyImg.type());
  78. //exclude the point which intersect into this bounding box but is not in this group.
  79. int label = i;
  80. for (int row = 0; row < rectROI.rows; row++)
  81. {
  82. for (int col = 0; col < rectROI.cols; col++)
  83. {
  84. int v = rectROI.at<int>(row, col);
  85. if (v == label)
  86. {
  87. imageROI.at<uchar>(row, col) = 255;//set the value to 255,so we won't consider other pixel when we find segment and feature in this ROI.
  88. }
  89. }
  90. }
  91. COTSParticleList roiParts;
  92. if (GetOneParticleFromROI(rectMax.x, rectMax.y, rectMax.width, rectMax.height, imageROI.data, roiParts))
  93. {
  94. if (roiParts.size() > 0)
  95. {
  96. COTSParticlePtr roiPart = roiParts[0];//we will find only one part in the roi.
  97. roiPart->SetXRayPos(CPoint(center_x, center_y));
  98. CRect r = CRect(x, y, x + w, y + h);
  99. roiPart->SetParticleRect(r);
  100. roiPart->SetActualArea(actualArea);
  101. roiPart->SetPixelArea(area);
  102. listParticleOut.push_back(roiPart);
  103. }
  104. }
  105. }
  106. }
  107. int nTagId;
  108. for (auto pParticle : listParticleOut)
  109. {
  110. COTSFeaturePtr pFeature = pParticle->GetFeature();
  111. COTSSegmentsList listSegment = pFeature->GetSegmentsList();
  112. long nPixelNum = 0;
  113. long nPixelAll = 0;
  114. int nStartS = 0;
  115. int nHeightS = 0;
  116. int nLengthS = 0;
  117. for (auto pSegment : listSegment)
  118. {
  119. // get particle average gray
  120. nStartS = pSegment->GetStart();
  121. nHeightS = pSegment->GetHeight();
  122. nLengthS = pSegment->GetLength();
  123. nPixelNum += (long)nLengthS;
  124. for (unsigned int i = 0; i < nLengthS; i++)
  125. {
  126. long nValueTemp = (long)*(pSrcImg + nHeightS * nWidthImg + nStartS + i);
  127. nPixelAll += nValueTemp;
  128. }
  129. }
  130. BYTE nAveGray = (BYTE)(nPixelAll / nPixelNum);
  131. pParticle->SetAveGray(nAveGray);
  132. auto fieldOTSRect = m_pFieldData->GetOTSRect();
  133. CPoint leftTop = fieldOTSRect.GetTopLeft();
  134. CRect rectInSinglefld = pParticle->GetParticleRect();
  135. CPoint OTSLeftTop = CPoint(leftTop.x + rectInSinglefld.left * a_pixelSize, leftTop.y - rectInSinglefld.top * a_pixelSize);
  136. CPoint OTSRightBottom = CPoint(leftTop.x + rectInSinglefld.right * a_pixelSize, leftTop.y - rectInSinglefld.bottom * a_pixelSize);
  137. COTSRect recInOTSCord = COTSRect(OTSLeftTop, OTSRightBottom);
  138. pParticle->SetOTSRect(recInOTSCord);
  139. }
  140. m_pFieldData->SetParticleList(listParticleOut);
  141. }
  142. delete[]pTempImg;
  143. return TRUE;
  144. }
  145. BOOL COTSImageProcess::GetParticlesBySpecialGrayRange(CBSEImgPtr a_pBSEImg, CIntRangePtr a_grayRange,CDoubleRangePtr a_diameterRange,double a_pixelSize, COTSFieldDataPtr m_pFieldData)
  146. {
  147. ASSERT(m_pFieldData);
  148. ASSERT(a_pBSEImg);
  149. ASSERT(a_grayRange);
  150. int nWidthImg = a_pBSEImg->GetWidth();
  151. int nHeightImg = a_pBSEImg->GetHeight();
  152. m_pFieldData->Width = nWidthImg;
  153. m_pFieldData->Height = nHeightImg;
  154. long nImgSize = nWidthImg * nHeightImg;
  155. BYTE* pSrcImg = a_pBSEImg->GetImageDataPointer();
  156. BYTE* pTempImg = new BYTE[nImgSize];
  157. CRect r = CRect(0, 0, nWidthImg, nHeightImg);
  158. CBSEImgPtr imgNoBGBinary = CBSEImgPtr(new CBSEImg(r));
  159. long nNumParticle = 0;
  160. GetSpecialGrayRangeImage(a_pBSEImg, a_grayRange, imgNoBGBinary, nNumParticle);
  161. BYTE* pPixel = imgNoBGBinary->GetImageDataPointer();
  162. if (nNumParticle == 0)
  163. {
  164. COTSParticleList listParticleEmpty;
  165. listParticleEmpty.clear();
  166. m_pFieldData->SetParticleList(listParticleEmpty);
  167. }
  168. else
  169. {
  170. // get the area image
  171. Mat cvcopyImg = Mat(nHeightImg, nWidthImg, CV_8UC1, pPixel);
  172. Mat labels = Mat::zeros(cvcopyImg.size(), CV_32S);
  173. Mat stats, centroids;
  174. int number = connectedComponentsWithStats(cvcopyImg, labels, stats, centroids, 8, CV_32S);
  175. double rStart = a_diameterRange->GetStart() / 2.0;
  176. double rEnd = a_diameterRange->GetEnd() / 2.0;
  177. double areaStart = rStart * rStart * 3.14159;
  178. double areaEnd = rEnd * rEnd * 3.14159;
  179. COTSParticleList listParticleOut;
  180. for (size_t i = 1; i < number; i++)
  181. {
  182. int center_x = centroids.at<double>(i, 0);
  183. int center_y = centroids.at<double>(i, 1);
  184. //앤近긋움
  185. int x = stats.at<int>(i, CC_STAT_LEFT);
  186. int y = stats.at<int>(i, CC_STAT_TOP);
  187. int w = stats.at<int>(i, CC_STAT_WIDTH);
  188. int h = stats.at<int>(i, CC_STAT_HEIGHT);
  189. int area = stats.at<int>(i, CC_STAT_AREA);
  190. double actualArea = area * a_pixelSize * a_pixelSize;
  191. if (actualArea >= areaStart && actualArea < areaEnd)
  192. {
  193. Rect rectMax = Rect(x, y, w, h);
  194. Mat rectROI = labels(rectMax).clone();
  195. Mat imageROI = Mat::zeros(rectMax.size(), cvcopyImg.type());
  196. //exclude the point which intersect into this bounding box but is not in this group.
  197. int label = i;
  198. for (int row = 0; row < rectROI.rows; row++)
  199. {
  200. for (int col = 0; col < rectROI.cols; col++)
  201. {
  202. int v = rectROI.at<int>(row, col);
  203. if (v == label)
  204. {
  205. imageROI.at<uchar>(row, col) = 255;
  206. }
  207. }
  208. }
  209. COTSParticleList roiParts;
  210. if (!GetOneParticleFromROI(rectMax.x, rectMax.y, rectMax.width, rectMax.height, imageROI.data, roiParts))
  211. {
  212. continue;
  213. }
  214. if (roiParts.size() > 0)
  215. {
  216. COTSParticlePtr roiPart = roiParts[0];
  217. roiPart->SetXRayPos(CPoint(center_x, center_y));
  218. CRect r = CRect(x, y, x + w, y + h);
  219. roiPart->SetParticleRect(r);
  220. roiPart->SetActualArea(actualArea);
  221. roiPart->SetPixelArea(area);
  222. listParticleOut.push_back(roiPart);
  223. }
  224. }
  225. }
  226. // form a image only have particles on
  227. //COTSSegmentsList listImage;
  228. for (auto pParticle : listParticleOut)
  229. {
  230. int area = pParticle->GetActualArea();
  231. double pActualArea = area ;
  232. COTSFeaturePtr pFeature = pParticle->GetFeature();
  233. COTSSegmentsList listSegment = pFeature->GetSegmentsList();
  234. long nPixelNum = 0;
  235. long nPixelAll = 0;
  236. int nStartS = 0;
  237. int nHeightS = 0;
  238. int nLengthS = 0;
  239. for (auto pSegment : listSegment)
  240. {
  241. // get particle average gray
  242. nStartS = pSegment->GetStart();
  243. nHeightS = pSegment->GetHeight();
  244. nLengthS = pSegment->GetLength();
  245. nPixelNum += (long)nLengthS;
  246. for (unsigned int i = 0; i < nLengthS; i++)
  247. {
  248. long nValueTemp = (long)*(pSrcImg + nHeightS * nWidthImg + nStartS + i);
  249. nPixelAll += nValueTemp;
  250. }
  251. }
  252. BYTE nAveGray = (BYTE)(nPixelAll / nPixelNum);
  253. pParticle->SetAveGray(nAveGray);
  254. auto fieldOTSRect = m_pFieldData->GetOTSRect();
  255. CPoint leftTop = fieldOTSRect.GetTopLeft();
  256. CRect rectInSinglefld = pParticle->GetParticleRect();
  257. CPoint OTSLeftTop = CPoint(leftTop.x + rectInSinglefld.left * a_pixelSize, leftTop.y - rectInSinglefld.top * a_pixelSize);
  258. CPoint OTSRightBottom = CPoint(leftTop.x + rectInSinglefld.right * a_pixelSize, leftTop.y - rectInSinglefld.bottom * a_pixelSize);
  259. COTSRect recInOTSCord = COTSRect(OTSLeftTop, OTSRightBottom);
  260. pParticle->SetOTSRect(recInOTSCord);
  261. }
  262. m_pFieldData->SetParticleList(listParticleOut);
  263. }
  264. delete[]pTempImg;
  265. return TRUE;
  266. }
  267. CIntRangePtr COTSImageProcess::CalBackground(CBSEImgPtr m_pBSEImg)
  268. {
  269. auto ranges = CalcuGrayLevelRange(m_pBSEImg);
  270. return ranges[0];
  271. }
  272. std::vector<CIntRangePtr> COTSImageProcess::CalcuGrayLevelRange(CBSEImgPtr m_pBSEImg)
  273. {
  274. CIntRangePtr pBackground = CIntRangePtr(new CIntRange());
  275. WORD originChartData[MAXBYTE];
  276. WORD firstSmoothChart[MAXBYTE];
  277. WORD secondSmooth[MAXBYTE];
  278. //1. get chart data
  279. m_pBSEImg->SetChartData();
  280. memcpy(originChartData, m_pBSEImg->GetBSEChart(), sizeof(WORD) * MAXBYTE);
  281. originChartData[0] = 0;
  282. originChartData[254] = 0;
  283. linearSmooth5(originChartData, firstSmoothChart, MAXBYTE);
  284. linearSmooth5(firstSmoothChart, secondSmooth, MAXBYTE);
  285. //2. get down edge
  286. int nLengthEdge = MAXBYTE + 2;
  287. WORD n_aBSEChart[MAXBYTE + 2];
  288. memset(n_aBSEChart, 0, sizeof(WORD) * nLengthEdge);
  289. std::map<long, std::vector<int>> peakMap;// hold all the peaks in this spectrum which are sorted by there area.
  290. std::vector<int> currentUpSeries;
  291. std::vector<int> currentPeakSeries;
  292. // make sure the wave begin with up edge and end with down edge
  293. n_aBSEChart[0] = 0;
  294. n_aBSEChart[nLengthEdge - 1] = 0;
  295. memcpy(&n_aBSEChart[1], &secondSmooth, sizeof(WORD) * MAXBYTE);
  296. int nLengthCom = MAXBYTE + 1;
  297. // up edge
  298. for (int i = 0; i < nLengthCom; i++)
  299. {
  300. if (n_aBSEChart[i] <= n_aBSEChart[i + 1])//this is a upward edge
  301. {
  302. if (currentPeakSeries.size() > 0)
  303. {
  304. int seriesSize = currentPeakSeries.size();
  305. long area = 0;
  306. for (int i = 0; i < seriesSize; i++)
  307. {
  308. area = area + n_aBSEChart[currentPeakSeries[i]];
  309. }
  310. peakMap[area] = currentPeakSeries;
  311. currentPeakSeries.clear();
  312. }
  313. currentUpSeries.push_back(i + 1);// save all the continuous up edge
  314. }
  315. else//this is a downward edge
  316. {
  317. // encounter a downward edge means upward edge series end,
  318. if (currentUpSeries.size() > 0)
  319. {
  320. currentPeakSeries = currentUpSeries;
  321. currentUpSeries.clear();
  322. }
  323. currentPeakSeries.push_back(i + 1);
  324. }
  325. }
  326. if (currentPeakSeries.size() > 0)
  327. {
  328. int seriesSize = currentPeakSeries.size();
  329. long area = 0;
  330. for (int i = 0; i < seriesSize; i++)
  331. {
  332. area = area + n_aBSEChart[currentPeakSeries[i]];
  333. }
  334. peakMap[area] = currentPeakSeries;
  335. currentPeakSeries.clear();
  336. }
  337. std::vector<CIntRangePtr> ranges;
  338. std::map<long, std::vector<int>>::reverse_iterator it;
  339. for (it=peakMap.rbegin();it!=peakMap.rend();it++)
  340. {
  341. CIntRangePtr pRange = CIntRangePtr(new CIntRange());
  342. pRange->SetStart(it->second[0]);
  343. pRange->SetEnd(it->second[it->second.size()-1]);
  344. ranges.push_back(pRange);
  345. }
  346. return ranges;
  347. }
  348. void COTSImageProcess::GetSpecialGrayRangeImage(CBSEImgPtr a_pImgIn, CIntRangePtr a_SpecialGrayRange, CBSEImgPtr a_pBinImgOut, long& foundedPixelNum)
  349. {
  350. // the background pixel will be 0,and the other part will be 255.
  351. ASSERT(a_pImgIn);
  352. int nWidthImg = a_pImgIn->GetWidth();
  353. int nHeightImg = a_pImgIn->GetHeight();
  354. long nImgSize = nWidthImg * nHeightImg;
  355. BYTE* pTempImg = new BYTE[nImgSize];
  356. BYTE* pSrcImg = a_pImgIn->GetImageDataPointer();
  357. BYTE* pPixel = new byte[nImgSize];
  358. long nBGStart;
  359. long nBGEnd;
  360. long nNumParticle = 0;
  361. nBGStart = a_SpecialGrayRange->GetStart();
  362. nBGEnd = a_SpecialGrayRange->GetEnd();
  363. // delete background
  364. for (unsigned int i = 0; i < nImgSize; i++)
  365. {
  366. if (pSrcImg[i] >= nBGStart && pSrcImg[i] <= nBGEnd)
  367. {
  368. pPixel[i] = 255;
  369. nNumParticle++;
  370. }
  371. else
  372. {
  373. pPixel[i] = 0;
  374. }
  375. }
  376. BErode3(pPixel, pTempImg, 5, nHeightImg, nWidthImg);
  377. BDilate3(pTempImg, pPixel, 5, nHeightImg, nWidthImg);
  378. a_pBinImgOut->SetImageData(pPixel, nWidthImg, nHeightImg);
  379. foundedPixelNum = nNumParticle;
  380. delete[] pTempImg;
  381. return;
  382. }
  383. void COTSImageProcess::RemoveBackGround(CBSEImgPtr a_pImgIn, COTSImageProcessParamPtr a_pImageProcessParam, CBSEImgPtr a_pBinImgOut,long& foundedPixelNum)
  384. {
  385. // the background pixel will be 0,and the other part will be 255.
  386. ASSERT(a_pImgIn);
  387. ASSERT(a_pImageProcessParam);
  388. int nWidthImg = a_pImgIn->GetWidth();
  389. int nHeightImg = a_pImgIn->GetHeight();
  390. long nImgSize = nWidthImg * nHeightImg;
  391. BYTE* pTempImg = new BYTE[nImgSize];
  392. BYTE* pSrcImg = a_pImgIn->GetImageDataPointer();
  393. BYTE* pPixel= new BYTE[nImgSize];
  394. Mat srcImgMat = GetMatDataFromBseImg(a_pImgIn);
  395. Mat rstMat;
  396. long nBGStart;
  397. long nBGEnd;
  398. long nPartStart;
  399. long nPartEnd;
  400. long nNumParticle = 0;
  401. if (a_pImageProcessParam->GetBGRemoveType() == OTS_BGREMOVE_TYPE::MANUAL)
  402. {
  403. nBGStart = a_pImageProcessParam->GetBGGray().GetStart();
  404. nBGEnd = a_pImageProcessParam->GetBGGray().GetEnd();
  405. nPartStart = a_pImageProcessParam->GetParticleGray().GetStart();
  406. nPartEnd = a_pImageProcessParam->GetParticleGray().GetEnd();
  407. CVRemoveBG(srcImgMat, rstMat, nBGStart, nBGEnd,nNumParticle);
  408. //RemoveBG_old(srcImgMat, rstMat, nBGStart, nBGEnd, nNumParticle);
  409. pPixel = rstMat.data;
  410. }
  411. else
  412. {
  413. auto range = CalBackground(a_pImgIn);
  414. nBGStart = range->GetStart();
  415. nBGEnd = range->GetEnd();
  416. switch (a_pImageProcessParam->GetAutoBGRemoveType())
  417. {
  418. case OTS_AUTOBGREMOVE_TYPE::DOWNWARD:
  419. //RemoveBG_old(srcImgMat, rstMat, 0, nBGEnd, nNumParticle);
  420. CVRemoveBG(srcImgMat, rstMat, 0, nBGEnd, nNumParticle);
  421. break;
  422. case OTS_AUTOBGREMOVE_TYPE::UPWARD:
  423. //RemoveBG_old(srcImgMat, rstMat, nBGStart, 255, nNumParticle);
  424. CVRemoveBG(srcImgMat, rstMat, nBGStart, 255, nNumParticle);
  425. break;
  426. case OTS_AUTOBGREMOVE_TYPE::MIDDLE:
  427. //RemoveBG_old(srcImgMat, rstMat, nBGStart, nBGEnd, nNumParticle);
  428. CVRemoveBG(srcImgMat, rstMat, nBGStart, nBGEnd, nNumParticle);
  429. break;
  430. default:
  431. break;
  432. }
  433. }
  434. a_pBinImgOut->SetImageData(pPixel,nWidthImg,nHeightImg);
  435. foundedPixelNum = nNumParticle;
  436. delete[] pTempImg;
  437. return ;
  438. }
  439. BOOL COTSImageProcess::GetParticles(long left, long top, long a_nWidth, long a_nHeight, const BYTE* a_pPixel, COTSParticleList& a_listParticles)
  440. {
  441. ASSERT(a_pPixel);
  442. if (!a_pPixel)
  443. {
  444. return FALSE;
  445. }
  446. //a_listParticles.clear();
  447. COTSParticleList findedParts;
  448. COTSSegmentsList listSegment;
  449. listSegment.clear();
  450. //1. get segment line by line
  451. if (!GetSegmentList(left, top, a_nWidth, a_nHeight, a_pPixel, listSegment))
  452. {
  453. return FALSE;
  454. }
  455. if ((int)listSegment.size() == 0)
  456. {
  457. return FALSE;
  458. }
  459. //2. save the temp feature
  460. COTSFeatureList listFeature;
  461. listFeature.clear();
  462. if (!GetFeatureList(listSegment, listFeature))//get every feature for all the particle,the complete feature.
  463. {
  464. return FALSE;
  465. }
  466. if ((int)listFeature.size() == 0)
  467. {
  468. return FALSE;
  469. }
  470. /*COTSParticleList listParticles;
  471. listParticles.clear();*/
  472. if (!ChangeFeaturelist(listFeature, findedParts))
  473. {
  474. return FALSE;
  475. }
  476. for (auto f : findedParts)
  477. {
  478. a_listParticles.push_back(f);
  479. }
  480. return TRUE;
  481. }
  482. BOOL COTSImageProcess::GetOneParticleFromROI(long left, long top, long a_nWidth, long a_nHeight, const BYTE* a_pPixel, COTSParticleList& a_listParticles)
  483. {
  484. ASSERT(a_pPixel);
  485. if (!a_pPixel)
  486. {
  487. return FALSE;
  488. }
  489. //a_listParticles.clear();
  490. COTSParticleList findedParts;
  491. COTSSegmentsList listSegment;
  492. listSegment.clear();
  493. //1. get segment line by line
  494. if (!GetSegmentList(left, top, a_nWidth, a_nHeight, a_pPixel, listSegment))
  495. {
  496. return FALSE;
  497. }
  498. if ((int)listSegment.size() == 0)
  499. {
  500. return FALSE;
  501. }
  502. //2. save the temp feature
  503. COTSFeatureList listFeature;
  504. listFeature.clear();
  505. COTSFeaturePtr fea = COTSFeaturePtr(new COTSFeature());
  506. fea->SetSegmentsList(listSegment);
  507. listFeature.push_back(fea);
  508. if ((int)listFeature.size() == 0)
  509. {
  510. return FALSE;
  511. }
  512. if (!ChangeFeaturelist(listFeature, findedParts))
  513. {
  514. return FALSE;
  515. }
  516. for (auto f : findedParts)
  517. {
  518. a_listParticles.push_back(f);
  519. }
  520. return TRUE;
  521. }
  522. BOOL COTSImageProcess::GetSegmentList(long left, long top, long a_nWidth, long a_nHeight, const BYTE* a_pPixel, COTSSegmentsList& a_listSegments)
  523. {
  524. ASSERT(a_pPixel);
  525. long nImgSize = a_nWidth * a_nHeight;
  526. a_listSegments.clear();
  527. //1. get segment line by line
  528. long nLine, nm, nn;
  529. long nStart = 0, nLength = 0;
  530. for (nLine = 0; nLine < a_nHeight; nLine++)
  531. {
  532. for (nm = 0; nm < a_nWidth; nm += (nLength + 1))
  533. {
  534. nLength = 0;
  535. // get start
  536. if (*(a_pPixel + nLine * a_nWidth + nm) != 0)
  537. {
  538. nStart = nm;
  539. nLength++;
  540. //get length
  541. for (nn = nm + 1; nn < a_nWidth; nn++)
  542. {
  543. // check if segment is over, break
  544. if (nLength != 0)
  545. {
  546. if (*(a_pPixel + nLine * a_nWidth + nn) == 0)
  547. break;
  548. }
  549. if (*(a_pPixel + nLine * a_nWidth + nn) != 0)
  550. {
  551. nLength++;
  552. }
  553. }
  554. // generate segment
  555. COTSSegmentPtr pSegment = COTSSegmentPtr(new COTSSegment(nLine + top, nStart + left, nLength));
  556. a_listSegments.push_back(pSegment);
  557. }
  558. else
  559. {
  560. continue;
  561. }
  562. }
  563. }
  564. if ((int)a_listSegments.size() == 0)
  565. {
  566. //LogErrorTrace(__FILE__, __LINE__, _T("no particle is found."));
  567. return FALSE;
  568. }
  569. return TRUE;
  570. }
  571. BOOL COTSImageProcess::GetFeatureList(COTSSegmentsList& a_listSegments, COTSFeatureList& a_listFeatures)
  572. {
  573. COTSSegmentsList listSegmentNew;
  574. std::map<long, COTSSegmentsList > mapOneLineSegments;
  575. for each (auto s in a_listSegments)
  576. {
  577. mapOneLineSegments[s->GetHeight()].push_back(s);//sorting all the segments base on the line number.
  578. }
  579. std::map<long, COTSSegmentsList >::iterator lineItr = mapOneLineSegments.begin();//find the highest line
  580. while (lineItr != mapOneLineSegments.end())
  581. {
  582. for (auto s = lineItr->second.begin(); s < lineItr->second.end(); )//find one segment of this line.
  583. {
  584. COTSSegmentPtr bottomSeg = *s;
  585. listSegmentNew.clear();
  586. listSegmentNew.push_back(*s);
  587. s = lineItr->second.erase(s);
  588. std::map<long, COTSSegmentsList >::iterator tempItr = lineItr;
  589. tempItr++;
  590. for (; tempItr != mapOneLineSegments.end(); tempItr++)//find all other lines of segments
  591. {
  592. if (tempItr->first - bottomSeg->GetHeight() > 1)
  593. {
  594. break;
  595. }
  596. for (auto nextLineSegment = tempItr->second.begin(); nextLineSegment < tempItr->second.end();)//find next line's all segments
  597. {
  598. if (((*nextLineSegment)->GetStart() - (bottomSeg->GetStart() + bottomSeg->GetLength())) > 1)
  599. {
  600. break;
  601. }
  602. if (bottomSeg->UpDownConection(**nextLineSegment))
  603. {
  604. listSegmentNew.push_back(*nextLineSegment);
  605. bottomSeg = *nextLineSegment;
  606. nextLineSegment = tempItr->second.erase(nextLineSegment);
  607. break;
  608. }
  609. if (tempItr->second.size() > 0)
  610. {
  611. nextLineSegment++;
  612. }
  613. else
  614. {
  615. break;
  616. }
  617. }
  618. }
  619. COTSFeaturePtr pFeature = COTSFeaturePtr(new COTSFeature());
  620. pFeature->SetSegmentsList(listSegmentNew);
  621. //check if this new feature is connected with other found feature.
  622. COTSSegmentPtr topSeg = listSegmentNew[0];//find the toppest segment of this new feature.
  623. COTSSegmentPtr bottomSegment = listSegmentNew[listSegmentNew.size() - 1];//find the lowest segment of this new feature.
  624. bool haveMerged = false;
  625. for each (auto f in a_listFeatures)
  626. {
  627. for (auto seg : f->GetSegmentsList())
  628. {
  629. if (bottomSegment->UpDownConection(*seg) || topSeg->UpDownConection(*seg))
  630. {
  631. COTSSegmentsList segs = f->GetSegmentsList();
  632. for (auto s : listSegmentNew)
  633. {
  634. segs.push_back(s);
  635. }
  636. f->SetSegmentsList(segs);
  637. haveMerged = true;
  638. break;
  639. }
  640. }
  641. if (haveMerged)
  642. {
  643. break;
  644. }
  645. }
  646. if (!haveMerged)
  647. {
  648. a_listFeatures.push_back(pFeature);
  649. }
  650. if (lineItr->second.size() == 0)
  651. {
  652. break;
  653. }
  654. }
  655. lineItr++;
  656. }
  657. return true;
  658. }
  659. BOOL COTSImageProcess::ChangeFeaturelist(COTSFeatureList& a_listFeatures, COTSParticleList& a_listParticle)
  660. {
  661. for (auto pFeature : a_listFeatures)
  662. {
  663. COTSParticlePtr pParticle = COTSParticlePtr(new COTSParticle());
  664. pParticle->SetFeature(pFeature);
  665. a_listParticle.push_back(pParticle);
  666. }
  667. if ((int)a_listParticle.size() == 0)
  668. {
  669. return FALSE;
  670. }
  671. return TRUE;
  672. }
  673. BOOL COTSImageProcess::CalcuParticleImagePropertes(COTSParticlePtr a_pOTSPart, double a_PixelSize)
  674. {
  675. //--------- convert this particle data to image data,construct an image only with this particle.------
  676. const int nExpand_Size = 3;
  677. const int nWhiteColor = 0;
  678. const int nThickness = 1;
  679. // lineType Type of the line
  680. const int nLineType = 8;
  681. // get rectangle of the particle
  682. CRect rect = a_pOTSPart->GetParticleRect();
  683. if (a_pOTSPart->GetActualArea() < 30 * a_PixelSize)// the particle is too small that openCV can't calculate a width value of it. Then we take the upright rect of the particle as it's minArea rect.
  684. {
  685. double w = 0, h = 0;
  686. w = (double)rect.Width()*a_PixelSize;
  687. h = (double)rect.Height()*a_PixelSize;
  688. a_pOTSPart->SetDMax(MAX(w, h));
  689. a_pOTSPart->SetDMin(MIN(w, h));
  690. a_pOTSPart->SetDMean((w + h) / 2);
  691. a_pOTSPart->SetFeretDiameter((w + h) / 2);
  692. a_pOTSPart->SetDElong(MAX(w, h));
  693. a_pOTSPart->SetPerimeter((w+h)*2);
  694. a_pOTSPart->SetDPerp(MIN(w, h));
  695. a_pOTSPart->SetDInscr(MIN(w, h));
  696. a_pOTSPart->SetMinWidth(w);
  697. a_pOTSPart->SetMinHeight(h);
  698. return true;
  699. }
  700. // calculate the particle image data size, expand 3 pixel at the edge
  701. Mat particleImage = Mat::zeros(rect.Height() + nExpand_Size , rect.Width() + nExpand_Size , CV_8U);
  702. // get the segment list
  703. COTSSegmentsList listSegment = a_pOTSPart->GetFeature()->GetSegmentsList();
  704. for (auto pSegment : listSegment)
  705. {
  706. int nStart = pSegment->GetStart() - rect.left + nExpand_Size;
  707. int nEnd = pSegment->GetStart() + pSegment->GetLength() - rect.left - 1 + nExpand_Size;
  708. int nHeight = pSegment->GetHeight() - rect.top + nExpand_Size;
  709. line(particleImage, Point(nStart, nHeight), Point(nEnd, nHeight), Scalar(nBlackColor), nThickness, nLineType);
  710. }
  711. //--------abstract the contour of the particle.
  712. vector<vector<Point>>contours;
  713. findContours(particleImage, contours, CV_RETR_EXTERNAL, CV_CHAIN_APPROX_NONE);
  714. if (contours.size()==0)// the particle is too odd that openCV can't find a contour of it. Then we take the upright rect of the particle as it's minArea rect.
  715. {
  716. double w = 0, h = 0;
  717. w = (double)rect.Width()*a_PixelSize;
  718. h = (double)rect.Height()*a_PixelSize;
  719. a_pOTSPart->SetDMax(MAX(w, h));
  720. a_pOTSPart->SetDMin(MIN(w, h));
  721. a_pOTSPart->SetDMean((w + h) / 2);
  722. a_pOTSPart->SetFeretDiameter((w + h) / 2);
  723. a_pOTSPart->SetDElong(MAX(w, h));
  724. a_pOTSPart->SetPerimeter((w + h) * 2);
  725. a_pOTSPart->SetDPerp(MIN(w, h));
  726. a_pOTSPart->SetDInscr(MIN(w, h));
  727. a_pOTSPart->SetMinWidth(w);
  728. a_pOTSPart->SetMinHeight(h);
  729. return true;
  730. }
  731. int imaxcontour = 0, imax = 0;
  732. for (unsigned int i = 0; i < contours.size(); i++) {
  733. int itmp = contourArea(contours[i]);
  734. if (imaxcontour < itmp) {
  735. imax = i;
  736. imaxcontour = itmp;
  737. }
  738. }
  739. vector<Point > listEdge = contours[imax];
  740. vector<vector<Point>>Outcontours;
  741. Outcontours.push_back(listEdge);
  742. //---------calculate the minimium rectangle
  743. auto rRect = cv::minAreaRect(listEdge);
  744. Point2f p[4];
  745. rRect.points(p);
  746. int D_MIN =getDistance(p[0], p[1]);
  747. int D_MinRecLen = 0;//minareaRect's length(the longger side).
  748. for (int j = 0; j <= 2; j++)
  749. {
  750. //line(cvContourImg, p[j], p[(j + 1) % 4], Scalar(100, 100, 0), 2);
  751. int d = getDistance(p[j], p[j + 1]);
  752. if (d < D_MIN)
  753. {
  754. D_MIN = d;
  755. }
  756. if (d > D_MinRecLen)
  757. {
  758. D_MinRecLen = d;
  759. }
  760. }
  761. a_pOTSPart->SetDMin(D_MIN*a_PixelSize);
  762. a_pOTSPart->SetMinWidth(rRect.size.width);
  763. a_pOTSPart->SetMinHeight(rRect.size.height);
  764. //----------calculate the perimeter
  765. double d = arcLength(listEdge, true);
  766. a_pOTSPart->SetPerimeter(d*a_PixelSize);
  767. //-----------calculate the Max Diameter. Find the min enclosing circle first ,then find the two farthest circle connected point.
  768. Point2f center; float radius;
  769. minEnclosingCircle(listEdge, center, radius);
  770. //circle(cvContourImg, center, radius, Scalar(100), 2);
  771. std::vector <Point> outContour = listEdge;
  772. std::vector <Point> rst;
  773. for (unsigned int k = 0; k < outContour.size(); k++)
  774. {
  775. Point p = outContour[k];
  776. double d = sqrt(pow((p.x - center.x), 2) + pow((p.y - center.y), 2));
  777. if (fabs(d - radius) < 0.01)
  778. {
  779. rst.push_back(p);
  780. }
  781. }
  782. double D_MAX = 0;
  783. Point lineDmax[2];
  784. for (unsigned int m = 0; m < rst.size(); m++)
  785. {
  786. Point p = rst[m];
  787. for (unsigned int n = m + 1; n < rst.size(); n++)
  788. {
  789. Point p1 = rst[n];
  790. double d = sqrt(powf((p.x - p1.x), 2) + powf((p.y - p1.y), 2));
  791. if (d > D_MAX)
  792. {
  793. D_MAX = d;
  794. lineDmax[0] = p;
  795. lineDmax[1] = p1;
  796. }
  797. }
  798. }
  799. a_pOTSPart->SetDMax(D_MAX*a_PixelSize);
  800. //--------calculate the D_PERP property using the D_MAX's two endpoints.
  801. std::vector<Point> curve1;
  802. std::vector<Point> curve2;
  803. for (unsigned int i = 0; i < outContour.size(); i++)
  804. {
  805. Point pt = outContour[i];
  806. bool start = false;
  807. int clockwise = Side(lineDmax[0], lineDmax[1], pt);// devide these points into two group ,separate into the two sides.
  808. if (clockwise > 0)
  809. {
  810. curve1.push_back(pt);
  811. }
  812. else
  813. {
  814. curve2.push_back(pt);
  815. }
  816. }
  817. double d_perp1 = 0, d_perp2 = 0;
  818. for (unsigned int i = 0; i < curve1.size(); i++)
  819. {
  820. double d = getDist_P2L(curve1[i], lineDmax[0], lineDmax[1]);
  821. if (d > d_perp1)
  822. {
  823. d_perp1 = d;
  824. }
  825. }
  826. for (unsigned int i = 0; i < curve2.size(); i++)
  827. {
  828. double d = getDist_P2L(curve2[i], lineDmax[0], lineDmax[1]);
  829. if (d > d_perp2)
  830. {
  831. d_perp2 = d;
  832. }
  833. }
  834. a_pOTSPart->SetDPerp((d_perp1 + d_perp2)*a_PixelSize);
  835. //----------find the diameter of max inscribed circle
  836. int r;
  837. Point inscribeCirclecenter;
  838. FindInnerCircleInContour(outContour, inscribeCirclecenter, r);
  839. //--------------------------------------------------------calculate the xraypos !
  840. CRect rec = a_pOTSPart->GetParticleRect();
  841. a_pOTSPart->SetXRayPos(CPoint(inscribeCirclecenter.x - nExpand_Size + rec.left - 1, inscribeCirclecenter.y - nExpand_Size + rec.top - 1));
  842. a_pOTSPart->SetDInscr(r * 2 * a_PixelSize);
  843. //---------------calculate the image other caracater: length/width realArea/minRectangeArea etc. we can use these propertes to do forward process.
  844. double minRectArea = D_MIN * D_MinRecLen*a_PixelSize*a_PixelSize;//離鬼棍쌈앤近충생
  845. double fillRatio = a_pOTSPart->GetActualArea() / minRectArea;//茄셥충생宅離鬼棍쌈앤近충생궐,that's the fill rate.
  846. double lengthWidthRatio;
  847. lengthWidthRatio = (double)D_MinRecLen / D_MIN;//낀욱궐
  848. //decide if this shape is a strip shape :if the lenthWidthRatio>2 then it is. if the lengthWidthRatio<2 and the areaRatio<0.5 then it is.
  849. bool isStripShape = false;
  850. double curveLength = 0;
  851. double D_MEAN=0;
  852. Moments mu;
  853. mu = moments(listEdge, false);
  854. int nx = mu.m10 / mu.m00;
  855. int ny = mu.m01 / mu.m00;
  856. //circle(cvcopyImg, Point(nx, ny), 1, (255), 1);
  857. Point ptCenter = Point((int)nx, (int)ny);
  858. if (pointPolygonTest(listEdge, ptCenter, false) != 1)// the center point doesn't contain in the contour, we think it as curve shape.
  859. {
  860. isStripShape = true;
  861. }
  862. /*if (lengthWidthRatio >= 2 )// in PartA software this is true,but IncA because of the GB definition the everage feret diameter is always the mean value of all the chord.
  863. {
  864. isStripShape = true;
  865. }*/
  866. if (fillRatio <= 0.4)// only when the fill rate is very low,we think it as a curve shape,then we choose the mean width as the feret diameter.
  867. {
  868. isStripShape = true;
  869. }
  870. if (isStripShape)
  871. {
  872. curveLength = a_pOTSPart->GetPerimeter()/2 - a_pOTSPart->GetDInscr()/2;// thinking this particle as a strip rectangle.the width is the max inscribe circle diameter/2.
  873. if (curveLength < D_MAX)
  874. {
  875. curveLength = D_MAX;
  876. }
  877. if (curveLength < MIN_DOUBLE_VALUE || a_pOTSPart->GetActualArea()<MIN_DOUBLE_VALUE)
  878. {
  879. D_MEAN = 0;
  880. }
  881. else
  882. {
  883. D_MEAN = a_pOTSPart->GetActualArea() / curveLength;
  884. }
  885. a_pOTSPart->SetDMean(D_MEAN*a_PixelSize);
  886. a_pOTSPart->SetFeretDiameter(D_MEAN*a_PixelSize);
  887. a_pOTSPart->SetDElong (curveLength*a_PixelSize);
  888. }
  889. else//it's a ball shape particle
  890. {
  891. curveLength = D_MAX;
  892. double ftd = 0, maxD = 0, minD = 0, dratio = 0;
  893. GetParticleAverageChord(outContour, a_PixelSize, ftd);
  894. a_pOTSPart->SetDMean(ftd);
  895. a_pOTSPart->SetFeretDiameter(ftd);
  896. a_pOTSPart->SetDElong(curveLength*a_PixelSize);
  897. }
  898. if (a_pOTSPart->GetAspectRatio() > 2)
  899. {
  900. double angle;
  901. if (rRect.size.width > rRect.size.height) // w > h
  902. {
  903. angle = abs(rRect.angle);
  904. }
  905. else
  906. {
  907. angle = 90.0 + abs(rRect.angle);
  908. }
  909. a_pOTSPart->SetOrientation(angle);
  910. }
  911. else
  912. {
  913. a_pOTSPart->SetOrientation(0);
  914. }
  915. return true;
  916. }
  917. void COTSImageProcess::ImshowImage(CBSEImgPtr img)
  918. {
  919. BYTE* data = img->GetImageDataPointer();
  920. //Mat cvImg;
  921. cv::Size s;
  922. s.width = img->GetImageSize().cx;
  923. s.height = img->GetImageSize().cy;
  924. Mat cvImg=Mat::zeros(s, CV_8U);
  925. cvImg.data = data;
  926. cv::imshow("dd", cvImg);
  927. cv::waitKey();
  928. }
  929. void COTSImageProcess::ImshowChartData(CBSEImgPtr img)
  930. {
  931. img->SetChartData();
  932. WORD* data = img->GetBSEChart();
  933. //Mat cvImg;
  934. cv::Size s;
  935. s.width = 255;
  936. s.height = 100;
  937. Mat cvImg = Mat::zeros(s, CV_8U);
  938. //cvImg.data = data;
  939. WORD nBSEChart[MAXBYTE];
  940. //1. get chart data
  941. linearSmooth5(data, nBSEChart, MAXBYTE);
  942. for (int i=1;i<255;i++)
  943. {
  944. line(cvImg, Point(i, 100-nBSEChart[i]), Point(i+1, 100-nBSEChart[i+1]), Scalar(nBlackColor), 1, 8);
  945. }
  946. cv::imshow("chart", cvImg);
  947. cv::waitKey();
  948. }
  949. BOOL COTSImageProcess::MergeBigBoundaryParticles(COTSFieldDataList allFields,double pixelSize,int scanFieldSize, CSize ResolutionSize, COTSParticleList& mergedParts)
  950. {
  951. class BorderPart
  952. {
  953. typedef std::shared_ptr<BorderPart> CBorderPartPtr;
  954. BorderPart(COTSParticlePtr p)
  955. {
  956. myPart = p;
  957. headerParticle = NULL;
  958. }
  959. public:
  960. COTSParticlePtr myPart;
  961. COTSParticle* headerParticle;//used to merge particles ,if this particle has been merged then this pointer will point to the first particle of these merged particles or else it's NULL.
  962. static std::vector<CBorderPartPtr> ConvertPartToBorderPart(COTSParticleList parts)
  963. {
  964. static std::map<COTSParticle*, CBorderPartPtr> allborderPart;
  965. std::vector<CBorderPartPtr> borderParts;
  966. for (auto p : parts)
  967. {
  968. if (allborderPart.find(p.get()) == allborderPart.end())
  969. {
  970. auto borderp = CBorderPartPtr(new BorderPart(p));
  971. borderParts.push_back(borderp);
  972. allborderPart[p.get()] = borderp;
  973. }
  974. else
  975. {
  976. borderParts.push_back(allborderPart[p.get()]);
  977. }
  978. }
  979. return borderParts;
  980. }
  981. };
  982. auto FldMgr = new CFieldMgr(scanFieldSize, ResolutionSize);
  983. std::map<COTSParticle*, COTSParticleList> mapMergeParticles;//hold up all the boundary connected particles. the pair's first is also the member of these particles.
  984. std::map<COTSParticle*, COTSSegmentsList> mapMergedSegments;//hold up all the segment's corresponding clone in the connected particles.
  985. for (auto centerfld : allFields)
  986. {
  987. // find neighbor field on the left.
  988. auto leftFld = FldMgr->FindNeighborField(allFields, centerfld, SORTING_DIRECTION::LEFT);
  989. if (leftFld != nullptr)
  990. {
  991. auto lParts = centerfld->GetLeftBorderedBigParticles();
  992. auto rParts = leftFld->GetRightBorderedBigParticles();
  993. auto leftParts = BorderPart::ConvertPartToBorderPart(lParts);
  994. auto rightParts = BorderPart::ConvertPartToBorderPart(rParts);
  995. for (auto leftp : leftParts)
  996. {
  997. for (auto rightp : rightParts)
  998. {
  999. if (leftp->myPart->IsConnected(rightp->myPart.get(), centerfld->Width, centerfld->Height, (int)SORTING_DIRECTION::LEFT))
  1000. {
  1001. if (leftp->headerParticle != NULL)
  1002. {
  1003. if (rightp->headerParticle == NULL)
  1004. {
  1005. rightp->headerParticle = leftp->headerParticle;
  1006. mapMergeParticles[leftp->headerParticle].push_back(rightp->myPart);
  1007. }
  1008. }
  1009. else
  1010. {
  1011. if (rightp->headerParticle != NULL)
  1012. {
  1013. leftp->headerParticle = rightp->myPart.get();
  1014. mapMergeParticles[rightp->myPart.get()].push_back(leftp->myPart);
  1015. }
  1016. else
  1017. {
  1018. leftp->headerParticle = leftp->myPart.get();
  1019. rightp->headerParticle = leftp->myPart.get();
  1020. mapMergeParticles[leftp->myPart.get()].push_back(rightp->myPart);
  1021. }
  1022. }
  1023. }
  1024. }
  1025. }
  1026. }
  1027. //find neighbor field on the upward
  1028. auto upFld = FldMgr->FindNeighborField(allFields, centerfld, SORTING_DIRECTION::UP);
  1029. if (upFld != nullptr)
  1030. {
  1031. auto topBorderParts = centerfld->GetTopBorderedBigParticles();
  1032. auto bottomBorderParts = upFld->GetBottomBorderedBigParticles();
  1033. auto upParts = BorderPart::ConvertPartToBorderPart(topBorderParts);
  1034. auto downParts = BorderPart::ConvertPartToBorderPart(bottomBorderParts);
  1035. for (auto upprt : upParts)
  1036. {
  1037. for (auto downprt : downParts)
  1038. {
  1039. if (upprt->myPart->IsConnected(downprt->myPart.get(), centerfld->Width, centerfld->Height, (int)SORTING_DIRECTION::UP))
  1040. {
  1041. if (upprt->headerParticle != NULL)
  1042. {
  1043. if (downprt->headerParticle == NULL)
  1044. {
  1045. downprt->headerParticle = upprt->headerParticle;
  1046. mapMergeParticles[upprt->headerParticle].push_back(downprt->myPart);
  1047. }
  1048. }
  1049. else
  1050. {
  1051. if (downprt->headerParticle != NULL)
  1052. {
  1053. upprt->headerParticle = downprt->headerParticle;
  1054. mapMergeParticles[downprt->myPart.get()].push_back(upprt->myPart);
  1055. }
  1056. else
  1057. {
  1058. upprt->headerParticle = upprt->myPart.get();
  1059. downprt->headerParticle = upprt->myPart.get();
  1060. mapMergeParticles[upprt->myPart.get()].push_back(downprt->myPart);
  1061. }
  1062. }
  1063. }
  1064. }
  1065. }
  1066. }
  1067. //find neighbor field on the downward.
  1068. auto downFld = FldMgr->FindNeighborField(allFields, centerfld,SORTING_DIRECTION::DOWN);
  1069. if (downFld != nullptr)
  1070. {
  1071. auto bottomParts = centerfld->GetBottomBorderedBigParticles();
  1072. auto topParts = downFld->GetTopBorderedBigParticles();
  1073. auto downParts = BorderPart::ConvertPartToBorderPart(bottomParts);
  1074. auto upParts= BorderPart::ConvertPartToBorderPart(topParts);
  1075. for (auto downprt : downParts)
  1076. {
  1077. for (auto upprt : upParts)
  1078. {
  1079. if (downprt->myPart->IsConnected(upprt->myPart.get(), centerfld->Width, centerfld->Height, (int)SORTING_DIRECTION::DOWN))
  1080. {
  1081. if (downprt->headerParticle != NULL)
  1082. {
  1083. if (upprt->headerParticle == NULL)
  1084. {
  1085. upprt->headerParticle = downprt->headerParticle;
  1086. mapMergeParticles[downprt->headerParticle].push_back(upprt->myPart);
  1087. }
  1088. }
  1089. else
  1090. {
  1091. if (upprt->headerParticle != NULL)
  1092. {
  1093. downprt->headerParticle = upprt->headerParticle;
  1094. mapMergeParticles[upprt->headerParticle].push_back(downprt->myPart);
  1095. }
  1096. else
  1097. {
  1098. downprt->headerParticle = downprt->myPart.get();
  1099. upprt->headerParticle = downprt->myPart.get();
  1100. mapMergeParticles[downprt->myPart.get()].push_back(upprt->myPart);
  1101. }
  1102. }
  1103. }
  1104. }
  1105. }
  1106. }
  1107. //find neighbor field on the right.
  1108. auto rightFld = FldMgr->FindNeighborField(allFields, centerfld, SORTING_DIRECTION::RIGHT);
  1109. if (rightFld != nullptr)
  1110. {
  1111. auto rParts = centerfld->GetRightBorderedBigParticles();
  1112. auto lParts = rightFld->GetLeftBorderedBigParticles();
  1113. auto rightParts = BorderPart::ConvertPartToBorderPart(rParts);
  1114. auto leftParts = BorderPart::ConvertPartToBorderPart(lParts);
  1115. for (auto rightprt : rightParts)
  1116. {
  1117. for (auto leftprt : leftParts)
  1118. {
  1119. if (rightprt->myPart->IsConnected(leftprt->myPart.get(), centerfld->Width, centerfld->Height, (int)SORTING_DIRECTION::RIGHT))
  1120. {
  1121. if (rightprt->headerParticle != NULL)
  1122. {
  1123. if (leftprt->headerParticle == NULL)
  1124. {
  1125. leftprt->headerParticle = rightprt->headerParticle;
  1126. mapMergeParticles[rightprt->headerParticle].push_back(leftprt->myPart);
  1127. }
  1128. }
  1129. else
  1130. {
  1131. if (leftprt->headerParticle != NULL)
  1132. {
  1133. rightprt->headerParticle = leftprt->headerParticle;
  1134. mapMergeParticles[leftprt->headerParticle].push_back(rightprt->myPart);
  1135. }
  1136. else
  1137. {
  1138. rightprt->headerParticle = rightprt->myPart.get();
  1139. leftprt->headerParticle = rightprt->myPart.get();
  1140. mapMergeParticles[rightprt->myPart.get()].push_back(leftprt->myPart);
  1141. }
  1142. }
  1143. }
  1144. }
  1145. }
  1146. }
  1147. }
  1148. /*for (auto particle : mapMergeParticles)
  1149. {
  1150. }*/
  1151. static int partTagId;
  1152. for (auto pair : mapMergeParticles)
  1153. {
  1154. struct EleAreaPercentage
  1155. {
  1156. EleAreaPercentage(double p, CElementChemistryPtr e)
  1157. {
  1158. areaPercentage = p;
  1159. eleData = e;
  1160. }
  1161. double areaPercentage;
  1162. CElementChemistryPtr eleData;
  1163. };
  1164. auto newPart = COTSParticlePtr(new COTSParticle());
  1165. COTSSegmentsList newSegs;
  1166. auto p = pair.first;
  1167. newPart->SetSEMPos(p->GetSEMPos());
  1168. //firstly,we sum up all the merged particles's area and get the represent string.
  1169. std::string partsStr = std::to_string(p->GetFieldId()) + ":" + std::to_string(p->GetAnalysisId());
  1170. double allPartArea = p->GetActualArea();//Get the first particle's area.
  1171. for (auto other : pair.second)// Get the total area of all these particles for the use of ele calcu.
  1172. {
  1173. partsStr += "," + std::to_string(other->GetFieldId()) + ":" + std::to_string(other->GetAnalysisId());//Get the subparticles string such as "1:1,2:1" etc.
  1174. allPartArea += other->GetActualArea();//Get other particle's area
  1175. }
  1176. // calculate all the new segment's position.
  1177. std::vector <COTSParticle*> allSubParts;
  1178. allSubParts.push_back(p);
  1179. for (auto other : pair.second)// Get the total area of all these particles for the use of ele calcu.
  1180. {
  1181. allSubParts.push_back(other.get());
  1182. }
  1183. for (auto subp : allSubParts)
  1184. {
  1185. int fid = subp->GetFieldId();
  1186. CPoint myFldPos;
  1187. for (auto f : allFields)//find this particle's filed.
  1188. {
  1189. if (f->GetId() == fid)
  1190. {
  1191. myFldPos = f->GetPosition();
  1192. }
  1193. }
  1194. int fldWidth = allFields[0]->Width;
  1195. int fldHeight = allFields[0]->Height;
  1196. CPoint fldLeftUpPos = CPoint(myFldPos.x + fldWidth / 2 , myFldPos.y + fldHeight / 2 );
  1197. for (auto s : subp->GetFeature()->GetSegmentsList())
  1198. {
  1199. COTSSegmentPtr newseg = COTSSegmentPtr(new COTSSegment());
  1200. newseg->SetStart(s->GetStart() + fldLeftUpPos.x);
  1201. newseg->SetHeight((0 - s->GetHeight()) + fldLeftUpPos.y);//the coordinate system of segment in a field is different with the OTS coordinate system.OTS system's y axis is upward positive ,yet the field is downward positive.
  1202. newseg->SetLength(s->GetLength());
  1203. newSegs.push_back(newseg);
  1204. }
  1205. }
  1206. COTSFeaturePtr newFeature = COTSFeaturePtr(new COTSFeature());
  1207. newFeature->SetSegmentsList(newSegs);
  1208. newPart->SetFeature(newFeature);
  1209. newPart->CalCoverRectFromSegment();
  1210. //second, we get all the element data and their area percentage .
  1211. std::map<std::string, std::vector<EleAreaPercentage>> mapEleData;
  1212. CPosXrayPtr pXray1 = p->GetXrayInfo();
  1213. if (pXray1 != nullptr)
  1214. {
  1215. for (auto ele : pXray1->GetElementQuantifyData())
  1216. {
  1217. mapEleData[ele->GetName().GetBuffer()].push_back(EleAreaPercentage(p->GetActualArea() / allPartArea, ele));
  1218. }
  1219. }
  1220. for (auto other : pair.second)
  1221. {
  1222. auto otherXray = other->GetXrayInfo();
  1223. if (otherXray != nullptr)
  1224. {
  1225. for (auto eledata : otherXray->GetElementQuantifyData())
  1226. {
  1227. mapEleData[eledata->GetName().GetBuffer()].push_back(EleAreaPercentage(other->GetActualArea() / allPartArea, eledata));
  1228. }
  1229. }
  1230. }
  1231. // third,we calculate all the element's new percentage data and get a new element chemistry list.
  1232. CElementChemistriesList newCheList;
  1233. for (auto eledata : mapEleData)
  1234. {
  1235. CElementChemistryPtr newEleche = CElementChemistryPtr(new CElementChemistry());
  1236. newEleche->SetName(CString(eledata.first.c_str()));
  1237. double newPercentage = 0;
  1238. for (auto d : eledata.second)
  1239. {
  1240. newPercentage += d.areaPercentage * d.eleData->GetPercentage();
  1241. }
  1242. newEleche->SetPercentage(newPercentage);
  1243. newCheList.push_back(newEleche);
  1244. }
  1245. CPosXrayPtr xray(new CPosXray());
  1246. xray->SetElementQuantifyData(newCheList);
  1247. newPart->SetXrayInfo(xray);
  1248. newPart->SetConnectedParticlesSequentialString(partsStr);
  1249. newPart->SetActualArea(allPartArea);
  1250. partTagId++;
  1251. newPart->SetParticleId(partTagId);
  1252. newPart->SetAnalysisId(partTagId);
  1253. std::string name = p->GetClassifyName();
  1254. newPart->SetClassifyName(name);
  1255. newPart->SetColor(p->GetColor());
  1256. mergedParts.push_back(newPart);
  1257. }
  1258. return true;
  1259. }
  1260. }