/*------------------------------------------------------------------------------*
 * File Name: pro_image_utils.c	 												*
 * Creation: 2005-07-20															*
 * Purpose: OriginPro image utility functions									*
 * Copyright (c) OriginLab Corp. 2010											*
 * All Rights Reserved															*
 * 																				*
 * Modification Log:															*
 * Derek 03/25/2010 ADD_A_FUNCTION_FOR_FINDING_OBJECT_IN_BINARY_MATRIX          *
 *------------------------------------------------------------------------------*/
 
// Right-click on the line below and select 'Open "Origin.h"' to open the Origin.h
// system header file.
#include <Origin.h>
#include <ocim.h>			/// Fisher 03/24/10	New image processing module
#include <ocStats.h>
#include <pro_image_utils.h>


// Derek 03/25/2010 ADD_A_FUNCTION_FOR_FINDING_OBJECT_IN_BINARY_MATRIX 
int find_objects_in_binary_matrix(const matrix<uint>& im, vector<uint>& ov, int minObjPixels, int maxObjPixels)
{
    int nRet;
	int nRows, nCols, nSize;
	
	nRows = im.GetNumRows();
	nCols = im.GetNumCols();
	nSize = nRows * nCols;
		

	// Label the connected components of the binary matrix
	vector<uint>	vLabel;
	vLabel.SetSize(nSize);
	int nRegNum = ocim_connected_component_label(im, nRows, nCols, vLabel);	
	if(nRegNum < 0)
	{
		return -1;
	}

	// Calculate the areas of the labeled regions
	vector<int> vAreas;
    vAreas.SetSize(nRegNum);
	nRet = ocim_object_feature_area(vLabel, nSize, nRegNum, vAreas);
	if(nRet < 0)
	{
		return -1;
	}
	
	// Count the occurrences of regions with area within range [minObjPixels, maxObjPixels]
    double minArea, maxArea;
	vAreas.GetMinMax(minArea, maxArea);
	if(maxObjPixels == -1) maxObjPixels = (int)maxArea;

	int nSz = maxObjPixels - minObjPixels + 1;
    ov.SetSize(nSz);
	
	nRet = ocmath_frequency_count_i(vAreas, (UINT)nRegNum, minObjPixels, maxObjPixels, ov);
    if (nRet < 0)
    {
		return -1;
	}
    
	double nObjNum;
	ov.Sum(nObjNum); // the total number of objects found

	return ((int)nObjNum);	
}
//end ADD_A_FUNCTION_FOR_FINDING_OBJECT_IN_BINARY_MATRIX 



int particle_analysis_instance(const matrix<double>& mSrcImage, matrix<uint>& mSegmLabels, matrix<uint>& mObjectLabels)
{

    int ii, jj;

	UINT nRows, nCols, nSize;	
	nRows = mSrcImage.GetNumRows();
	nCols = mSrcImage.GetNumCols();
	nSize = nRows * nCols;
    const double nMaxPixValue = mSrcImage.GetMax();    


	// Initial parameters for the GMM segmentation algorithm	

	UINT nClassNum = 3; // This parameter must be specified manually

	GMM1DModel initModel;
	initModel.nCompNum = nClassNum;
	initModel.pMu      = (double *)malloc(nClassNum*sizeof(double));
	initModel.pSigma   = (double *)malloc(nClassNum*sizeof(double));
	initModel.pWeight  = (double *)malloc(nClassNum*sizeof(double));
	for(ii = 0; ii < nClassNum; ii ++)
	{
		initModel.pMu[ii]     = nMaxPixValue * ii / (double)nClassNum;
		initModel.pSigma[ii]  = 1.0;
		initModel.pWeight[ii] = 1.0 / nClassNum;
	}

	
	// Image segmentation
	ocim_segmentation_by_gmm1d(mSrcImage, nSize, nClassNum, mSegmLabels, &initModel, 2);

	
	// Object feature parameters (must be specified manaully)
	ObjectFeatures ObjFeatures;	    
	ObjFeatures.bIsAreaPresent = true;
	ObjFeatures.nMinArea       = 3;
	ObjFeatures.nMaxArea       = 10000;
	ObjFeatures.bIsOrientationPresent = false;
	ObjFeatures.fMinOrientAng       = 0;
	ObjFeatures.fMaxOrientAng       = 10;	 
    
	// Identify objects by features
	int nRegNum = ocim_identify_object_by_features(mSegmLabels, nRows, nCols, &ObjFeatures, mObjectLabels);
       
	
	return nRegNum;
}


void particle_analysis_instance_run()
{
	
	MatrixPage SrcImgPage = Project.MatrixPages(0); 
	if(!SrcImgPage)
	{
		return;
	}	
	MatrixLayer SrcImgLayer(SrcImgPage.GetName());    	
    MatrixObject mo = SrcImgLayer.MatrixObjects();
	
    matrix<double> mSrcImg(mo);
	    
	UINT nRows, nCols, nSize;	
	nRows = mSrcImg.GetNumRows();
	nCols = mSrcImg.GetNumCols();
	nSize = nRows * nCols;
	    
	  
    matrix<uint> mSegmLabels;
    mSegmLabels.SetSize(nRows, nCols);
    
    matrix<uint> mObjectLabels;	
	mObjectLabels.SetSize(nRows, nCols);
	
	UINT nRegNum = particle_analysis_instance(mSrcImg, mSegmLabels, mObjectLabels);    
    
	//printf("Identified Particle Number = %d\n\n", nRegNum);
        
        
    MatrixPage SegmLabelsPage;
    SegmLabelsPage.Create("SegmLabels");
    MatrixLayer SegmLabelsLayer(SegmLabelsPage.GetName());
    
    if( SegmLabelsLayer )
    {
       Matrix    ma(SegmLabelsLayer);
       ma = mSegmLabels;       
    }

   	
    MatrixPage ObjectLabelsPage;
    ObjectLabelsPage.Create("ObjectLabels");
    MatrixLayer ObjectLabelsLayer(ObjectLabelsPage.GetName());
   
    if( ObjectLabelsLayer )
    {
       Matrix    mb(ObjectLabelsLayer);
       mb = mObjectLabels;       
    }
    
    SegmLabelsLayer.SetViewImage(true);
    ObjectLabelsLayer.SetViewImage(true);
    
    	
	char str[200];
	sprintf(str, "Identified Particle Number = %d\n\n", nRegNum);
	MessageBox( GetWindow(), str, "Result");
}
