OTS/OTSCPP/OTSClassifyEngine/CurveCompareEngine/OTSClassifyOnCurveCompEng.cpp

#pragma once
#include "stdafx.h"
#include "OTSSTDLib.h"
#include "OTSClassifyOnCurveCompEng.h"
#include "OTSHelper.h"

using namespace std;
using namespace OTSTools;
namespace OTSClassifyEngine
{
	

	CClassifyOnCurveCompEng::CClassifyOnCurveCompEng()						// constructor
	{
	}

	CClassifyOnCurveCompEng::~CClassifyOnCurveCompEng()						// detractor
	{
	}
	int CClassifyOnCurveCompEng::GetCompareMineral(int iStartChannel, const vector<CString>& vecIgnoreElementNames, CPosXrayPtr posxray, CSTDLibPtr& oretype, double& dValue)
	{
		if (oretype->GetSTDItemCount() < 1)
		{
			return -1;
		}
		int id = 0;	// 存储最相似的那个mineralid
		double dcos = 0;
		dValue = -100000;	// 获取最大的相似值
		for (unsigned int i = 0; i < oretype->GetSTDItemCount(); i++)
		{
			dcos = GetCosValue(iStartChannel, vecIgnoreElementNames, posxray, oretype->GetSTDItem(i)->GetXrayData(), oretype->GetSTDItem (i)->GetChannelsNo());
			if (dcos > dValue)
			{
				dValue = dcos;
				id = oretype->GetSTDItem (i)->GetID();
			}
		}

		return id;
	}
	CStandardItemPtr CClassifyOnCurveCompEng::GetCompareSTD(int iStartChannel, const vector<CString>& vecIgnoreElementNames, CPosXrayPtr posxray, CSTDLibPtr& oretype, double& dValue)
	{
		CStandardItemPtr itm;	// 存储最相似的那个mineral
 		if (oretype->GetSTDItemCount() < 1)
		{
			return itm;
		}
		double dcos = 0;
		dValue = -100000;	// 获取最大的相似值
		for (unsigned int i = 0; i < oretype->GetSTDItemCount(); i++)
		{
			dcos = GetCosValue(iStartChannel, vecIgnoreElementNames, posxray, oretype->GetSTDItem(i)->GetXrayData(), oretype->GetSTDItem(i)->GetChannelsNo());
			if (dcos > dValue)
			{
				dValue = dcos;
				itm = oretype->GetSTDItem(i);
			}
		}

		return itm;
	}

	// 20190903:增加参数iStartChannel哪个位置开始比较, vecIgnoreChannel里面的通道不处理,
	double CClassifyOnCurveCompEng::GetCosValue(int iStartChannel, const vector<CString>& vecIgnoreElementNames, CPosXrayPtr posxray, DWORD* pXrayData1, int iDataLen1)
	{
		if (posxray->GetChannelsNo() != iDataLen1)
		{
			return 0;
		}

		vector<int> vecIgnoreChannel;
		GetElementChannel(vecIgnoreElementNames, vecIgnoreChannel);

		if (iStartChannel < 0 || iStartChannel >= iDataLen1)
		{
			return 1;
		}

		bool bignore = false;
		int j = 0;

		// 公式: (x1y1+x2y2+x3y3+...x2000y2000) / (sqrt(x1^2 + x2^2 + ...x2000^2) * sqrt(y1^2 + y2^2 + ...y2000^2))
		double dotpro = 0;
		double d1 = 0;
		double d2 = 0;
		for (int i = iStartChannel; i < iDataLen1; i++)
		{
			bignore = false;
			for (j = 0; j < (int)vecIgnoreChannel.size(); j++)
			{
				if (vecIgnoreChannel[j] == i)
				{
					bignore = true;
					break;
				}
			}

			if (!bignore)
			{
				dotpro = dotpro + posxray->GetXrayData()[i] * pXrayData1[i];
				d1 = d1 + posxray->GetXrayData()[i] * posxray->GetXrayData()[i];
				d2 = d2 + pXrayData1[i] * pXrayData1[i];
			}

		}
		d1 = sqrt(d1);
		d2 = sqrt(d2);
//            return (0 == d1 || 0 == d2) ? 0 : dotpro / (d1 * d2);
		// 算法改进, 加上距离权重
		if (0 == d1 || 0 == d2)
		{
			return 0;
		}
		double dresult = 0;
		dresult = (d1 < d2) ? d1 / d2 : d2 / d1;
		return 0.4 * dotpro / (d1 * d2) + 0.3 * dresult + 0.3 * GetStdEvp(iStartChannel, vecIgnoreChannel, posxray, pXrayData1, iDataLen1);
	}

	// 20190903:增加参数iStartChannel哪个位置开始比较, vecIgnoreChannel里面的通道不处理
	double CClassifyOnCurveCompEng::GetStdEvp(int iStartChannel, const vector<int>& vecIgnoreChannel, CPosXrayPtr posxray, DWORD* pXrayData1, int iDataLen1)
	{
		if (posxray->GetChannelsNo() != iDataLen1)
		{
			return 0;
		}
		if (iStartChannel < 0 || iStartChannel >= iDataLen1)
		{
			return 1;
		}

		bool bignore = false;
		int i = 0;
		int j = 0;
		int iSum1 = 0;
		int iSum2 = 0;
		double dAva1 = 0;
		double dAva2 = 0;
		double dSquare1 = 0;
		double dSquare2 = 0;
		vector<int> vecdata1;
		vector<int> vecdata2;

		vecdata1.clear();
		vecdata2.clear();

		// 求出能谱的平均值		
		for (i = iStartChannel; i < iDataLen1; i++)
		{
			bignore = false;
			for (j = 0; j < (int)vecIgnoreChannel.size(); j++)
			{
				if (vecIgnoreChannel[j] == i)
				{
					bignore = true;
					break;
				}
			}

			// 能谱差之和
			if (!bignore)
			{
				iSum1 = iSum1 + (int)pXrayData1[i];
				vecdata1.push_back((int)pXrayData1[i]);

				iSum2 = iSum2 + (int)posxray->GetXrayData()[i];
				vecdata2.push_back((int)posxray->GetXrayData()[i]);
			}
		}

		dAva1 = iSum1 / vecdata1.size();
		dAva2 = iSum2 / vecdata2.size();

		for (i = 0; i < (int)vecdata1.size(); i++)
		{
			dSquare1 = dSquare1 + (vecdata1[i] - dAva1) * (vecdata1[i] - dAva1);
			dSquare2 = dSquare2 + (vecdata2[i] - dAva1) * (vecdata2[i] - dAva1);
		}
		dSquare1 = dSquare1 / vecdata1.size();
		dSquare2 = dSquare2 / vecdata2.size();

		dSquare1 = sqrt(dSquare1);
		dSquare2 = sqrt(dSquare2);

		dSquare1 = dSquare1 * (double)iSum2 / iSum1;    // 归一化

		return (dSquare1 < dSquare2) ? (dSquare1 / dSquare2) : (dSquare2 / dSquare1);
	}

	void CClassifyOnCurveCompEng::GetElementChannel(const vector<CString>& vecstrElementNames, vector<int>& veciChannel)
	{
		veciChannel.clear();
		for (unsigned int i = 0; i < (int)vecstrElementNames.size(); i++)
		{
			GetElementChannel(vecstrElementNames[i], veciChannel);
		}
	}

	void CClassifyOnCurveCompEng::GetElementChannel(CString strElementName, vector<int>& veciChannel)
	{
		CElement* pElement = NULL;
		pElement = new CElement(strElementName);
		if (pElement)
		{
			AddChannel(pElement->GetEnergyValueK(), veciChannel);
			AddChannel(pElement->GetEnergyValueL(), veciChannel);
			AddChannel(pElement->GetEnergyValueM(), veciChannel);

			delete pElement;
			pElement = NULL;
		}
	}

	void CClassifyOnCurveCompEng::AddChannel(double dEnergy, vector<int>& veciChannel)
	{
		int i = 0;
		int iChannel = 0;
		if (dEnergy >= 0 && dEnergy < 20)
		{
			dEnergy = dEnergy * 100;
			iChannel = (int)dEnergy;
			for (i = 0; i < (int)veciChannel.size(); i++)
			{
				if (veciChannel[i] == iChannel)
				{
					break;
				}
			}
			if (i >= (int)veciChannel.size())
			{
				veciChannel.push_back(iChannel);
			}
		}
	}
}