opencv EM

 void myEM(int clusterNum, vector<Point>samplesPtVec, int cov_mat_type, vector<int> &LabelsVec,Mat &mean,vector<Mat> &varMatVec)
 {
	 int nsamplesNum = samplesPtVec.size();
	 Mat samples( nsamplesNum, 2, CV_32FC1 );

	 for( int i = 0; i < nsamplesNum; i++ )
	 {
		 samples.at<float>(i,0) = samplesPtVec[i].x;
		 samples.at<float>(i,1) = samplesPtVec[i].y;
	 }

	 Mat labels;
	 CvEM em_model;
	 CvEMParams params;
	 // initialize model parameters
	 params.covs      = NULL;
	 params.means     = NULL;
	 params.weights   = NULL;
	 params.probs     = NULL;
	 params.nclusters = clusterNum;
	 params.cov_mat_type  = cov_mat_type;//CvEM::COV_MAT_GENERIC;//COV_MAT_DIAGONAL COV_MAT_GENERIC  COV_MAT_SPHERICAL
	 params.start_step         = CvEM::START_AUTO_STEP;
	 params.term_crit.max_iter = 500;
	 params.term_crit.epsilon   = 0.01;
	 params.term_crit.type       = CV_TERMCRIT_ITER|CV_TERMCRIT_EPS;

	 // cluster the data
	 em_model.train(samples, Mat(), params, &labels);
	 mean = em_model.getMeans();

	 //cout<<" ++++++++++++++++++"<<endl;
	 //cout<<" mean = "<<mean<<endl;
	 //vector<Mat> varMatVec;
	 em_model.getCovs(varMatVec);
	 // 	for (int i=0;i<varMatVec.size();i++)
	 // 	{
	 // 		cout<<" var = "<<varMatVec[i]<<endl;
	 // 		Mat ss,uu,vt;
	 // 		cv::SVDecomp(varMatVec[i],ss,uu,vt,SVD::FULL_UV);  
	 // 		cout<<" ss = "<<ss<<endl;
	 // 		cout<<" uu = "<<uu<<endl;
	 // 		cout<<" vt  = "<<vt<<endl;
	 // 	}

	 LabelsVec.clear();
	 LabelsVec.resize(nsamplesNum);

	 for( int i = 0; i < nsamplesNum; i++ )
	 {
		 LabelsVec[i] = labels.at<int>(i);
	 }
 }


 void drawEM(Mat rgbImg, Mat mean, vector<Mat> varMatVec,Scalar color)
 {
	 int clusterNum = varMatVec.size();

	 for (int i=0;i<clusterNum;i++)
	 {
		 Point cPt = Point(mean.at<double>(i,0),mean.at<double>(i,1));
		 cv::circle(rgbImg,cPt,4,color/*colors[i]*/,2);

		 Mat ss,uu,vt;
		 cv::SVDecomp(varMatVec[i],ss,uu,vt,SVD::FULL_UV);  

		 Size sss = Size(sqrt(ss.at<double>(0,0))*1.5, sqrt(ss.at<double>(1,0))*1.5);
		 double angle = atan(uu.at<double>(1,0)/(uu.at<double>(0,0)+1/1000000))*180/CV_PI;

		 // 		if (sss.width > sss.height*1.2)
		 // 		{
		 // 			cv::ellipse(rgbImg, cPt,sss, angle, 0, 360, color, 3);
		 // 		}
		 // 		else
		 {
			 cv::ellipse(rgbImg, cPt,sss, angle, 0, 360, color, 2);//CV_AA
		 }
		 cout<<" sss = "<<sss<<" angle = "<<angle<<endl;
	 }
 }


void drawPtsEM(Mat rgbImg,vector<Point> ptVec,Scalar color)
{
	Mat mean;
	vector<Mat> varMatVec;
	vector<int>labelsVec;
	//COV_MAT_SPHERICAL  COV_MAT_DIAGONAL COV_MAT_GENERIC
	myEM(1, ptVec,CvEM::COV_MAT_GENERIC, labelsVec, mean, varMatVec);
	drawEM(rgbImg,mean,varMatVec,color);
	//cout<<" varMatVec = "<<varMatVec[0]<<endl;
}