OTSClassifyOnSpectrumCompEng.cpp 5.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234
  1. #pragma once
  2. #include "stdafx.h"
  3. #include "OTSSTDLib.h"
  4. #include "OTSClassifyOnSpectrumCompEng.h"
  5. #include "OTSHelper.h"
  6. using namespace std;
  7. using namespace OTSTools;
  8. namespace OTSClassifyEngine
  9. {
  10. CClassifyOnSpectrumCompEng::CClassifyOnSpectrumCompEng(int iStartChannel, CSTDLibPtr stdLib ) // constructor
  11. {
  12. m_iStartChannel = iStartChannel;
  13. m_stdLib = stdLib;
  14. }
  15. CClassifyOnSpectrumCompEng::~CClassifyOnSpectrumCompEng()// detractor
  16. {
  17. }
  18. int CClassifyOnSpectrumCompEng::GetMatchingSTDId( CPosXrayPtr posxray, double& dValue)
  19. {
  20. if (m_stdLib->GetSTDItemCount() < 1)
  21. {
  22. return -1;
  23. }
  24. int id = 0; // 存储最相似的那个mineralid
  25. double dcos = 0;
  26. dValue = -100000; // 获取最大的相似值
  27. for (unsigned int i = 0; i < m_stdLib->GetSTDItemCount(); i++)
  28. {
  29. dcos = GetCosValue( posxray, m_stdLib->GetSTDItem(i)->GetXraySpectrum(), m_stdLib->GetSTDItem (i)->GetXraySpectrum()->GetChannelsNum());
  30. if (dcos > dValue)
  31. {
  32. dValue = dcos;
  33. id = m_stdLib->GetSTDItem (i)->GetID();
  34. }
  35. }
  36. return id;
  37. }
  38. CSTDSpectrumItemPtr CClassifyOnSpectrumCompEng::GetMatchingSTD( CPosXrayPtr posxray, double& dValue)
  39. {
  40. CSTDSpectrumItemPtr itm; // 存储最相似的那个mineral
  41. int stdItemCount = m_stdLib->GetSTDItemCount();
  42. if (stdItemCount < 1)
  43. {
  44. return itm;
  45. }
  46. double dcos = 0;
  47. double dSimilarity = -1000;
  48. // 获取最大的相似值
  49. int idx = -1;
  50. for (unsigned int i = 0; i < stdItemCount; i++)
  51. {
  52. dcos = GetCosValue(posxray, m_stdLib->GetSTDItem(i)->GetXraySpectrum(), m_stdLib->GetSTDItem(i)->GetXraySpectrum()->GetChannelsNum());
  53. if (dcos > dSimilarity)
  54. {
  55. dSimilarity = dcos;
  56. idx = i;
  57. }
  58. }
  59. dValue = dSimilarity;
  60. itm = m_stdLib->GetSTDItem(idx);
  61. return itm;
  62. }
  63. CSTDSpectrumItemList CClassifyOnSpectrumCompEng::GetSTDOrderListBySpectrumMatching(CPosXrayPtr spectrum)
  64. {
  65. CSTDSpectrumItemList itms;
  66. std::map<double, CSTDSpectrumItemPtr> similaritymap;
  67. int stdItemCount = m_stdLib->GetSTDItemCount();
  68. if (stdItemCount < 1)
  69. {
  70. return itms;
  71. }
  72. double dcos = 0;
  73. double dSimilarity = -1000;
  74. // 获取最大的相似值
  75. int idx = -1;
  76. for (unsigned int i = 0; i < stdItemCount; i++)
  77. {
  78. dcos = GetCosValue(spectrum, m_stdLib->GetSTDItem(i)->GetXraySpectrum(), m_stdLib->GetSTDItem(i)->GetXraySpectrum()->GetChannelsNum());
  79. similaritymap[dcos] = m_stdLib->GetSTDItem(i);
  80. }
  81. auto itmIter = similaritymap.rbegin();
  82. while (itmIter!=similaritymap.rend())
  83. {
  84. itms.push_back(itmIter->second);
  85. itmIter++;
  86. }
  87. return itms;
  88. }
  89. // 20190903:增加参数iStartChannel哪个位置开始比较, vecIgnoreChannel里面的通道不处理,
  90. double CClassifyOnSpectrumCompEng::GetCosValue( CPosXrayPtr posXray,CPosXrayPtr posXray1, int iDataLen1)
  91. {
  92. if (posXray->GetChannelsNum() != iDataLen1)
  93. {
  94. return 0;
  95. }
  96. DWORD* pXrayData = posXray->GetXrayData();
  97. DWORD* pXrayData1 = posXray1->GetXrayData();
  98. // 公式: (x1y1+x2y2+x3y3+...x2000y2000) / (sqrt(x1^2 + x2^2 + ...x2000^2) * sqrt(y1^2 + y2^2 + ...y2000^2))
  99. double dotProduct = 0;
  100. double d1 = 0;
  101. double d2 = 0;
  102. for (int i = m_iStartChannel; i < iDataLen1; i++)
  103. {
  104. double r1 = pXrayData[i];
  105. double r2= pXrayData1[i];
  106. r1 *= r2;
  107. dotProduct = dotProduct + r1;
  108. }
  109. d1 = posXray->GetXrayDataVectorNorm();
  110. d2 = posXray1->GetXrayDataVectorNorm();
  111. return (0 == d1 || 0 == d2) ? 0 : dotProduct / (d1 * d2);
  112. }
  113. // 20190903:增加参数iStartChannel哪个位置开始比较, vecIgnoreChannel里面的通道不处理
  114. double CClassifyOnSpectrumCompEng::GetStdEvp( CPosXrayPtr posxray, DWORD* pXrayData1, int iDataLen1)
  115. {
  116. if (posxray->GetChannelsNum() != iDataLen1)
  117. {
  118. return 0;
  119. }
  120. DWORD* pXrayData = posxray->GetXrayData();
  121. //bool bignore = false;
  122. DWORD iSum1 = 0;
  123. DWORD iSum2 = 0;
  124. double dAva1 = 0;
  125. double dAva2 = 0;
  126. double dSquare1 = 0;
  127. double dSquare2 = 0;
  128. int num1= iDataLen1 - m_iStartChannel;
  129. // 求出能谱的平均值
  130. for (int i = m_iStartChannel; i < iDataLen1; i++)
  131. {
  132. // 能谱差之和
  133. iSum1 = iSum1 + pXrayData1[i];
  134. iSum2 = iSum2 + pXrayData[i];
  135. }
  136. dAva1 = iSum1 /num1;
  137. dAva2 = iSum2 / num1;
  138. for (int i = m_iStartChannel; i < iDataLen1; i++)
  139. {
  140. dSquare1 = dSquare1 + (pXrayData1[i] - dAva1) * (pXrayData1[i] - dAva1);
  141. dSquare2 = dSquare2 + (pXrayData[i] - dAva1) * (pXrayData[i] - dAva1);
  142. }
  143. dSquare1 = dSquare1 / num1;
  144. dSquare2 = dSquare2 / num1;
  145. dSquare1 = sqrt(dSquare1);
  146. dSquare2 = sqrt(dSquare2);
  147. dSquare1 = dSquare1 * dAva2 / dAva1; // 归一化
  148. return (dSquare1 < dSquare2) ? (dSquare1 / dSquare2) : (dSquare2 / dSquare1);
  149. }
  150. void CClassifyOnSpectrumCompEng::GetElementChannel(const vector<CString>& vecstrElementNames, vector<int>& veciChannel)
  151. {
  152. veciChannel.clear();
  153. for (unsigned int i = 0; i < (int)vecstrElementNames.size(); i++)
  154. {
  155. GetElementChannel(vecstrElementNames[i], veciChannel);
  156. }
  157. }
  158. void CClassifyOnSpectrumCompEng::GetElementChannel(CString strElementName, vector<int>& veciChannel)
  159. {
  160. CElement* pElement = NULL;
  161. pElement = new CElement(strElementName);
  162. if (pElement)
  163. {
  164. AddChannel(pElement->GetEnergyValueK(), veciChannel);
  165. AddChannel(pElement->GetEnergyValueL(), veciChannel);
  166. AddChannel(pElement->GetEnergyValueM(), veciChannel);
  167. delete pElement;
  168. pElement = NULL;
  169. }
  170. }
  171. void CClassifyOnSpectrumCompEng::AddChannel(double dEnergy, vector<int>& veciChannel)
  172. {
  173. int i = 0;
  174. int iChannel = 0;
  175. if (dEnergy >= 0 && dEnergy < 20)
  176. {
  177. dEnergy = dEnergy * 100;
  178. iChannel = (int)dEnergy;
  179. for (i = 0; i < (int)veciChannel.size(); i++)
  180. {
  181. if (veciChannel[i] == iChannel)
  182. {
  183. break;
  184. }
  185. }
  186. if (i >= (int)veciChannel.size())
  187. {
  188. veciChannel.push_back(iChannel);
  189. }
  190. }
  191. }
  192. }