SMS partial tracking bug fixes

3 files changed, 746 insertions, 559 deletions

diff --git a/sms.py b/sms.py
index 4eb53c0..1780998 100644
--- a/sms.py
+++ b/sms.py
@@ -37,13 +37,14 @@ class SMSPeakDetection(simpl.PeakDetection):
         self._analysis_params.iFrameRate = self.sampling_rate / self._hop_size
         self._analysis_params.iWindowType = pysms.SMS_WIN_HAMMING
         self._analysis_params.fHighestFreq = 20000
-        self._analysis_params.iMaxDelayFrames = 21#4
+        self._analysis_params.iMaxDelayFrames = 4
         self._analysis_params.analDelay = 0
         self._analysis_params.minGoodFrames = 1
         self._analysis_params.iCleanTracks = 0
         self._analysis_params.iFormat = pysms.SMS_FORMAT_HP
         self._analysis_params.nTracks = self._max_peaks
         self._analysis_params.maxPeaks = self._max_peaks
+        #self._analysis_params.nGuides = self._max_peaks
         pysms.sms_initAnalysis(self._analysis_params)
         self._peaks = pysms.SMS_SpectralPeaks(self.max_peaks)
         # By default, SMS will change the size of the frames being read depending on the
@@ -185,11 +186,10 @@ class SMSPartialTracking(simpl.PartialTracking):
         self._analysis_params = pysms.SMS_AnalParams()
         self._analysis_params.iSamplingRate = self.sampling_rate
         self._analysis_params.fHighestFreq = 20000
-        self._analysis_params.fLowestFundamental = 50
-        self._analysis_params.fDefaultFundamental = 100
-        self._analysis_params.iMaxDelayFrames = 3 # minimum frame delay with libsms
+        self._analysis_params.iMaxDelayFrames = 4 # minimum frame delay with libsms
         self._analysis_params.analDelay = 0
         self._analysis_params.minGoodFrames = 1
+        self._analysis_params.iCleanTracks = 0
         self._analysis_params.iFormat = pysms.SMS_FORMAT_HP
         self._analysis_params.nTracks = self.max_partials
         self._analysis_params.nGuides = self.max_partials
diff --git a/sms/peakContinuation.c b/sms/peakContinuation.c
index 7af4388..571772f 100644
--- a/sms/peakContinuation.c
+++ b/sms/peakContinuation.c
@@ -30,104 +30,100 @@
 #define GUIDE_ACTIVE 0
 
 #define MAX_CONT_CANDIDATES 5  /*!< maximum number of peak continuation
-                                  candidates */
+                                 candidates */
 /*! \brief function to get the next closest peak from a guide
  *
- * \param fGuideFreq		guide's frequency
- * \param pFFreqDistance	distance of last best peak from guide
- * \param pSpectralPeaks	array of peaks
- * \param pAnalParams	        analysis parameters
- * \param fFreqDev		        maximum deviation from guide
+ * \param fGuideFreq        guide's frequency
+ * \param pFFreqDistance    distance of last best peak from guide
+ * \param pSpectralPeaks    array of peaks
+ * \param pAnalParams           analysis parameters
+ * \param fFreqDev              maximum deviation from guide
  * \return peak number or -1 if nothing is good
  */
-int GetNextClosestPeak (sfloat fGuideFreq, sfloat *pFFreqDistance,
-                               SMS_Peak *pSpectralPeaks, SMS_AnalParams *pAnalParams,
-                               sfloat fFreqDev)
+int GetNextClosestPeak(sfloat fGuideFreq, sfloat *pFFreqDistance,
+                       SMS_Peak *pSpectralPeaks, SMS_AnalParams *pAnalParams,
+                       sfloat fFreqDev)
 {
-	int iInitialPeak = 
-		SMS_MAX_NPEAKS * fGuideFreq / (pAnalParams->iSamplingRate * .5),
-		iLowPeak, iHighPeak, iChosenPeak = -1;
-	sfloat fLowDistance, fHighDistance, fFreq;
-
-	if (pSpectralPeaks[iInitialPeak].fFreq <= 0)
-		iInitialPeak = 0;
-	  
-	/* find a low peak to start */
-//	if(fGuideFreq > 825 && fGuideFreq < 827)
-//	{
-//		printf("%d\t%f\t%f\n", SMS_MAX_NPEAKS, *pFFreqDistance,fGuideFreq - pSpectralPeaks[iInitialPeak].fFreq);
-//	}
-	fLowDistance = fGuideFreq - pSpectralPeaks[iInitialPeak].fFreq;
-	if (floor(fLowDistance) < floor(*pFFreqDistance))
-	{
-		while (floor(fLowDistance) <= floor(*pFFreqDistance) && 
-		       iInitialPeak > 0)
-		{
-			iInitialPeak--;
-			fLowDistance = fGuideFreq - pSpectralPeaks[iInitialPeak].fFreq;
-		}
-	}
-	else
-	{
-		while (floor(fLowDistance) >= floor(*pFFreqDistance) &&
-		       iInitialPeak < SMS_MAX_NPEAKS)
-		{
-			iInitialPeak++;	 
-//			if(fGuideFreq > 825 && fGuideFreq < 827)
-//			{
-//				printf("%d, %f\n", iInitialPeak, pSpectralPeaks[iInitialPeak].fFreq);
-//			}
-			if ((fFreq = pSpectralPeaks[iInitialPeak].fFreq) == 0)
-				return -1;
-
-			fLowDistance = fGuideFreq - fFreq;
-		}
-		iInitialPeak--;
-		fLowDistance = fGuideFreq - pSpectralPeaks[iInitialPeak].fFreq;
-	}
-
-	if (floor(fLowDistance) <= floor(*pFFreqDistance) || 
-	    fLowDistance > fFreqDev)
-		iLowPeak = -1;
-	else
-		iLowPeak = iInitialPeak;
-    
-	/* find a high peak to finish */
-	iHighPeak = iInitialPeak;
-	fHighDistance = fGuideFreq - pSpectralPeaks[iHighPeak].fFreq;
-	while (floor(fHighDistance) >= floor(-*pFFreqDistance) &&
-	       iHighPeak < SMS_MAX_NPEAKS)
-	{
-		iHighPeak++;	 
-		if ((fFreq = pSpectralPeaks[iHighPeak].fFreq) == 0)
-		{
-			iHighPeak = -1;
-			break;
-		}
-		fHighDistance = fGuideFreq - fFreq;
-	}
-	if (fHighDistance > 0 || fabs(fHighDistance) > fFreqDev ||
-	    floor(fabs(fHighDistance)) <= floor(*pFFreqDistance))
-		iHighPeak = -1;
-
-	/* chose between the two extrema */
-	if (iHighPeak >= 0 && iLowPeak >= 0)
-	{
-		if (fabs(fHighDistance) > fLowDistance)
-			iChosenPeak = iLowPeak;
-		else
-			iChosenPeak = iHighPeak;
-	}
-	else if (iHighPeak < 0 && iLowPeak >= 0)
-		iChosenPeak = iLowPeak;
-	else if (iHighPeak >= 0 && iLowPeak < 0)
-		iChosenPeak = iHighPeak;
-	else
-		return (-1);
-
-	*pFFreqDistance = fabs (fGuideFreq - pSpectralPeaks[iChosenPeak].fFreq);
-	return (iChosenPeak);
-}	
+    int iInitialPeak = SMS_MAX_NPEAKS * fGuideFreq / (pAnalParams->iSamplingRate * .5);
+    int iLowPeak, iHighPeak, iChosenPeak = -1;
+    sfloat fLowDistance, fHighDistance, fFreq;
+
+    if(pSpectralPeaks[iInitialPeak].fFreq <= 0)
+        iInitialPeak = 0;
+
+    if(iInitialPeak >= pAnalParams->maxPeaks)
+        iInitialPeak = 0;
+
+    fLowDistance = fGuideFreq - pSpectralPeaks[iInitialPeak].fFreq;
+    if(floor(fLowDistance) < floor(*pFFreqDistance))
+    {
+        while(floor(fLowDistance) <= floor(*pFFreqDistance) && 
+              iInitialPeak > 0)
+        {
+            iInitialPeak--;
+            fLowDistance = fGuideFreq - pSpectralPeaks[iInitialPeak].fFreq;
+        }
+    }
+    else
+    {
+        while(floor(fLowDistance) >= floor(*pFFreqDistance) &&
+                iInitialPeak < pAnalParams->maxPeaks)
+        {
+            iInitialPeak++; 
+            /* TODO: is this really the correct behaviour? Will ignore
+             * iInitialPeak values of maxPeaks */
+            if((fFreq = pSpectralPeaks[iInitialPeak].fFreq) == 0 || 
+               (iInitialPeak == pAnalParams->maxPeaks))
+               return -1;
+
+            fLowDistance = fGuideFreq - fFreq;
+        }
+        iInitialPeak--;
+        fLowDistance = fGuideFreq - pSpectralPeaks[iInitialPeak].fFreq;
+    }
+
+    if(floor(fLowDistance) <= floor(*pFFreqDistance) || 
+       fLowDistance > fFreqDev)
+        iLowPeak = -1;
+    else
+        iLowPeak = iInitialPeak;
+
+    /* find a high peak to finish */
+    iHighPeak = iInitialPeak;
+    fHighDistance = fGuideFreq - pSpectralPeaks[iHighPeak].fFreq;
+    while(floor(fHighDistance) >= floor(-*pFFreqDistance) &&
+          iHighPeak < pAnalParams->maxPeaks)
+    {
+        iHighPeak++;     
+        if((fFreq = pSpectralPeaks[iHighPeak].fFreq) == 0)
+        {
+            iHighPeak = -1;
+            break;
+        }
+        fHighDistance = fGuideFreq - fFreq;
+    }
+    if(fHighDistance > 0 || fabs(fHighDistance) > fFreqDev ||
+       floor(fabs(fHighDistance)) <= floor(*pFFreqDistance))
+        iHighPeak = -1;
+
+    /* chose between the two extrema */
+    if (iHighPeak >= 0 && iLowPeak >= 0)
+    {
+        if (fabs(fHighDistance) > fLowDistance)
+            iChosenPeak = iLowPeak;
+        else
+            iChosenPeak = iHighPeak;
+    }
+    else if (iHighPeak < 0 && iLowPeak >= 0)
+        iChosenPeak = iLowPeak;
+    else if (iHighPeak >= 0 && iLowPeak < 0)
+        iChosenPeak = iHighPeak;
+    else
+        return -1;
+
+    *pFFreqDistance = fabs (fGuideFreq - pSpectralPeaks[iChosenPeak].fFreq);
+    return iChosenPeak;
+}   
 
 /*! \brief choose the best candidate out of all
  *
@@ -136,378 +132,370 @@ int GetNextClosestPeak (sfloat fGuideFreq, sfloat *pFFreqDistance,
  * \param fFreqDev             maximum frequency deviation allowed
  * \return the peak number of the best candidate
  */
-static int ChooseBestCand (SMS_ContCandidate *pCandidate, int nCandidates, 
-                           sfloat fFreqDev)
+static int ChooseBestCand(SMS_ContCandidate *pCandidate, int nCandidates, sfloat fFreqDev)
 {
-	int i, iHighestCand, iClosestCand, iBestCand = 0;
-	sfloat fMaxMag, fClosestFreq;
-  
-	/* intial guess */
-	iClosestCand = 0;
-	fClosestFreq = pCandidate[iClosestCand].fFreqDev;
-	iHighestCand = 0;
-	fMaxMag = pCandidate[iHighestCand].fMagDev;
-  
-	/* get the best candidate */
-	for (i = 1; i < nCandidates; i++)
-	{
-		/* look for the one with highest magnitude */
-		if (pCandidate[i].fMagDev > fMaxMag)
-		{
-			fMaxMag = pCandidate[i].fMagDev;
-			iHighestCand = i;
-		}
-		/* look for the closest one to the guide */
-		if (pCandidate[i].fFreqDev < fClosestFreq)
-		{
-			fClosestFreq = pCandidate[i].fFreqDev;
-			iClosestCand = i;
-		}
-	}
-	iBestCand = iHighestCand;
-  
-	/* reconcile the two results */
-	if (iBestCand != iClosestCand &&
-	    fabs(pCandidate[iHighestCand].fFreqDev - fClosestFreq) > fFreqDev / 2)
-		iBestCand = iClosestCand; 
-  
-	return(pCandidate[iBestCand].iPeak);
+    int i, iHighestCand, iClosestCand, iBestCand = 0;
+    sfloat fMaxMag, fClosestFreq;
+
+    /* intial guess */
+    iClosestCand = 0;
+    fClosestFreq = pCandidate[iClosestCand].fFreqDev;
+    iHighestCand = 0;
+    fMaxMag = pCandidate[iHighestCand].fMagDev;
+
+    /* get the best candidate */
+    for (i = 1; i < nCandidates; i++)
+    {
+        /* look for the one with highest magnitude */
+        if (pCandidate[i].fMagDev > fMaxMag)
+        {
+            fMaxMag = pCandidate[i].fMagDev;
+            iHighestCand = i;
+        }
+        /* look for the closest one to the guide */
+        if (pCandidate[i].fFreqDev < fClosestFreq)
+        {
+            fClosestFreq = pCandidate[i].fFreqDev;
+            iClosestCand = i;
+        }
+    }
+    iBestCand = iHighestCand;
+
+    /* reconcile the two results */
+    if(iBestCand != iClosestCand &&
+       fabs(pCandidate[iHighestCand].fFreqDev - fClosestFreq) > fFreqDev / 2)
+        iBestCand = iClosestCand; 
+
+    return pCandidate[iBestCand].iPeak;
 }
 
 /*! \brief check for one guide that has choosen iBestPeak
  *
- * \param iBestPeak	choosen peak for a guide
- * \param pGuides		array of guides
- * \param nGuides		total number of guides
+ * \param iBestPeak choosen peak for a guide
+ * \param pGuides       array of guides
+ * \param nGuides       total number of guides
  * \return number of guide that chose the peak, or -1 if none
  */
-static int CheckForConflict (int iBestPeak, SMS_Guide *pGuides, int nGuides)
+static int CheckForConflict(int iBestPeak, SMS_Guide *pGuides, int nGuides)
 {
-	int iGuide;
-  
-	for (iGuide = 0; iGuide < nGuides; iGuide++)
-		if (pGuides[iGuide].iPeakChosen == iBestPeak)
-			return iGuide;
-   
-	return -1;
+    int iGuide;
+
+    for (iGuide = 0; iGuide < nGuides; iGuide++)
+        if (pGuides[iGuide].iPeakChosen == iBestPeak)
+            return iGuide;
+
+    return -1;
 }
 
 /*! \brief chose the best of the two guides for the conflicting peak 
  *
- * \param iConflictingGuide	conflicting guide number
- * \param iGuide			guide number
- * \param pGuides		        array of guides
- * \param pSpectralPeaks	array of peaks
+ * \param iConflictingGuide conflicting guide number
+ * \param iGuide            guide number
+ * \param pGuides               array of guides
+ * \param pSpectralPeaks    array of peaks
  * \return number of guide
  */
-static int BestGuide (int iConflictingGuide, int iGuide, SMS_Guide *pGuides,
-                      SMS_Peak *pSpectralPeaks)
+static int BestGuide(int iConflictingGuide, int iGuide, SMS_Guide *pGuides,
+                     SMS_Peak *pSpectralPeaks)
 {
-	int iConflictingPeak = pGuides[iConflictingGuide].iPeakChosen;
-	sfloat fGuideDistance = fabs (pSpectralPeaks[iConflictingPeak].fFreq -
-	                             pGuides[iGuide].fFreq);
-	sfloat fConfGuideDistance = fabs (pSpectralPeaks[iConflictingPeak].fFreq -
-	                                 pGuides[iConflictingGuide].fFreq);
-
-	if (fGuideDistance > fConfGuideDistance)
-		return (iConflictingGuide);
-	else
-		return (iGuide);
+    int iConflictingPeak = pGuides[iConflictingGuide].iPeakChosen;
+    sfloat fGuideDistance = fabs(pSpectralPeaks[iConflictingPeak].fFreq -
+                                 pGuides[iGuide].fFreq);
+    sfloat fConfGuideDistance = fabs(pSpectralPeaks[iConflictingPeak].fFreq -
+                                     pGuides[iConflictingGuide].fFreq);
+
+    if(fGuideDistance > fConfGuideDistance)
+        return iConflictingGuide;
+    else
+        return iGuide;
 }
 
 /*! \brief function to find the best continuation peak for a given guide
- * \param pGuides		guide attributes
- * \param iGuide		number of guide
- * \param pSpectralPeaks	peak values at the current frame
- * \param pAnalParams	        analysis parameters
+ * \param pGuides       guide attributes
+ * \param iGuide        number of guide
+ * \param pSpectralPeaks    peak values at the current frame
+ * \param pAnalParams           analysis parameters
  * \param fFreqDev                  frequency deviation allowed
  * \return the peak number
  */
-int GetBestPeak (SMS_Guide *pGuides, int iGuide, SMS_Peak *pSpectralPeaks,
-                        SMS_AnalParams *pAnalParams, sfloat fFreqDev)
+int GetBestPeak(SMS_Guide *pGuides, int iGuide, SMS_Peak *pSpectralPeaks,
+                SMS_AnalParams *pAnalParams, sfloat fFreqDev)
 {
-	int iCand = 0, iPeak, iBestPeak, iConflictingGuide, iWinnerGuide;
-	sfloat fGuideFreq = pGuides[iGuide].fFreq,
-		fGuideMag = pGuides[iGuide].fMag,
-		fMagDistance = 0;
-	sfloat fFreqDistance = -1;
-	SMS_ContCandidate pCandidate[MAX_CONT_CANDIDATES];
-
-	/* find all possible candidates */
-	while (iCand < MAX_CONT_CANDIDATES)
-	{
-		/* find the next best peak */
-		if ((iPeak = GetNextClosestPeak (fGuideFreq, &fFreqDistance,
-		                                 pSpectralPeaks, pAnalParams, 
-		                                 fFreqDev)) < 0)
-			break;
-	
-		/* if the peak's magnitude is not too small accept it as */
-		/* possible candidate        */
-		if ((fMagDistance = pSpectralPeaks[iPeak].fMag  - fGuideMag) > -20.0)
-		{
-			pCandidate[iCand].fFreqDev = fabs(fFreqDistance);
-			pCandidate[iCand].fMagDev = fMagDistance;
-			pCandidate[iCand].iPeak = iPeak;
-      	      
-			if(pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
-			   pAnalParams->iDebugMode == SMS_DBG_ALL)
-				fprintf (stdout, "candidate %d: freq %f mag %f\n", 
-				         iCand, pSpectralPeaks[iPeak].fFreq, 	
-				         pSpectralPeaks[iPeak].fMag);
-			iCand++;
-		}
-	}
-	/* get best candidate */
-	if (iCand < 1)
-		return (0);
-	else if (iCand == 1)
-		iBestPeak = pCandidate[0].iPeak;
-	else
-		iBestPeak = ChooseBestCand (pCandidate, iCand, 
-		                            pAnalParams->fFreqDeviation);
-      
-	if(pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
-	   pAnalParams->iDebugMode == SMS_DBG_ALL)
-		fprintf (stdout, "BestCandidate: freq %f\n",
-		         pSpectralPeaks[iBestPeak].fFreq);
-
-	/* if peak taken by another guide resolve conflict */
-	if ((iConflictingGuide = CheckForConflict (iBestPeak, pGuides, 
-	                                           pAnalParams->nGuides)) >= 0)
-	{
-		iWinnerGuide = BestGuide (iConflictingGuide, iGuide, pGuides, 
-		                          pSpectralPeaks);
-		if(pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
-		   pAnalParams->iDebugMode == SMS_DBG_ALL)
-			fprintf (stdout, 
-			         "Conflict: guide: %d (%f), and guide: %d (%f). best: %d\n", 
-			         iGuide, pGuides[iGuide].fFreq, 
-			         iConflictingGuide, pGuides[iConflictingGuide].fFreq, 	    
-			         iWinnerGuide);
-
-		if (iGuide == iWinnerGuide)				     
-		{
-			pGuides[iGuide].iPeakChosen = iBestPeak;
-			pGuides[iConflictingGuide].iPeakChosen = -1;
-		}
-	}
-	else
-		pGuides[iGuide].iPeakChosen = iBestPeak;
-
-	return (iBestPeak);
+    int iCand = 0, iPeak, iBestPeak, iConflictingGuide, iWinnerGuide;
+    sfloat fGuideFreq = pGuides[iGuide].fFreq,
+           fGuideMag = pGuides[iGuide].fMag,
+           fMagDistance = 0;
+    sfloat fFreqDistance = -1;
+    SMS_ContCandidate pCandidate[MAX_CONT_CANDIDATES];
+
+    /* find all possible candidates */
+    while (iCand < MAX_CONT_CANDIDATES)
+    {
+        /* find the next best peak */
+        if((iPeak = GetNextClosestPeak(fGuideFreq, &fFreqDistance,
+                                       pSpectralPeaks, pAnalParams, fFreqDev)) < 0)
+            break;
+
+        /* if the peak's magnitude is not too small accept it as */
+        /* possible candidate        */
+        if ((fMagDistance = pSpectralPeaks[iPeak].fMag  - fGuideMag) > -20.0)
+        {
+            pCandidate[iCand].fFreqDev = fabs(fFreqDistance);
+            pCandidate[iCand].fMagDev = fMagDistance;
+            pCandidate[iCand].iPeak = iPeak;
+
+            if(pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
+               pAnalParams->iDebugMode == SMS_DBG_ALL)
+                fprintf(stdout, "candidate %d: freq %f mag %f\n", 
+                        iCand, pSpectralPeaks[iPeak].fFreq,    
+                        pSpectralPeaks[iPeak].fMag);
+            iCand++;
+        }
+    }
+    /* get best candidate */
+    if(iCand < 1)
+        return 0;
+    else if (iCand == 1)
+        iBestPeak = pCandidate[0].iPeak;
+    else
+        iBestPeak = ChooseBestCand (pCandidate, iCand, 
+                pAnalParams->fFreqDeviation);
+
+    if(pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
+            pAnalParams->iDebugMode == SMS_DBG_ALL)
+        fprintf (stdout, "BestCandidate: freq %f\n",
+                pSpectralPeaks[iBestPeak].fFreq);
+
+    /* if peak taken by another guide resolve conflict */
+    if ((iConflictingGuide = CheckForConflict (iBestPeak, pGuides, 
+                    pAnalParams->nGuides)) >= 0)
+    {
+        iWinnerGuide = BestGuide (iConflictingGuide, iGuide, pGuides, 
+                pSpectralPeaks);
+        if(pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
+                pAnalParams->iDebugMode == SMS_DBG_ALL)
+            fprintf (stdout, 
+                    "Conflict: guide: %d (%f), and guide: %d (%f). best: %d\n", 
+                    iGuide, pGuides[iGuide].fFreq, 
+                    iConflictingGuide, pGuides[iConflictingGuide].fFreq,       
+                    iWinnerGuide);
+
+        if (iGuide == iWinnerGuide)                  
+        {
+            pGuides[iGuide].iPeakChosen = iBestPeak;
+            pGuides[iConflictingGuide].iPeakChosen = -1;
+        }
+    }
+    else
+        pGuides[iGuide].iPeakChosen = iBestPeak;
+
+    return iBestPeak;
 }
 
 /*! \brief function to get the next maximum (magnitude) peak
- * \param pSpectralPeaks	array of peaks
- * \param pFCurrentMax		last peak maximum
+ * \param pSpectralPeaks    array of peaks
+ * \param pFCurrentMax      last peak maximum
  * \return the number of the maximum peak
  */
 static int GetNextMax (SMS_Peak *pSpectralPeaks, sfloat *pFCurrentMax)
 {
-	sfloat fPeakMag;
-        sfloat fMaxMag = 0.;
-	int iPeak, iMaxPeak = -1;
-  
-	for (iPeak = 0; iPeak < SMS_MAX_NPEAKS; iPeak++)
-	{
-		fPeakMag = pSpectralPeaks[iPeak].fMag;
-    
-		if (fPeakMag == 0)
-			break;
-      
-		if (fPeakMag > fMaxMag && fPeakMag < *pFCurrentMax)
-		{
-			iMaxPeak = iPeak;
-			fMaxMag = fPeakMag;
-		}
-	}
-	*pFCurrentMax = fMaxMag;
-	return (iMaxPeak);
+    sfloat fPeakMag;
+    sfloat fMaxMag = 0.;
+    int iPeak, iMaxPeak = -1;
+
+    for (iPeak = 0; iPeak < SMS_MAX_NPEAKS; iPeak++)
+    {
+        fPeakMag = pSpectralPeaks[iPeak].fMag;
+
+        if (fPeakMag == 0)
+            break;
+
+        if (fPeakMag > fMaxMag && fPeakMag < *pFCurrentMax)
+        {
+            iMaxPeak = iPeak;
+            fMaxMag = fPeakMag;
+        }
+    }
+    *pFCurrentMax = fMaxMag;
+    return (iMaxPeak);
 }
 
 /*! \brief function to get a good starting peak for a track
  *
- * \param iGuide      		current guide
- * \param pGuides  		array of guides
- * \param nGuides			total number of guides
- * \param pSpectralPeaks	array of peaks
- * \param pFCurrentMax		current peak maximum
+ * \param iGuide            current guide
+ * \param pGuides       array of guides
+ * \param nGuides           total number of guides
+ * \param pSpectralPeaks    array of peaks
+ * \param pFCurrentMax      current peak maximum
  * \return \todo should this return something?
  */
-static int GetStartingPeak (int iGuide, SMS_Guide *pGuides, int nGuides,
-                            SMS_Peak *pSpectralPeaks, sfloat *pFCurrentMax)
+static int GetStartingPeak(int iGuide, SMS_Guide *pGuides, int nGuides,
+                           SMS_Peak *pSpectralPeaks, sfloat *pFCurrentMax)
 {
-	int iPeak = -1;
-	short peakNotFound = 1;
-  
-	while (peakNotFound == 1 && *pFCurrentMax > 0)
-	{
-                /* \todo I don't think this ever returns -1, but check */
-		if ((iPeak = GetNextMax (pSpectralPeaks, pFCurrentMax)) < 0)
-			return (-1);
-  
-		if (CheckForConflict (iPeak, pGuides, nGuides) < 0)
-		{
-			pGuides[iGuide].iPeakChosen = iPeak;
-			pGuides[iGuide].iStatus = GUIDE_BEG;
-			pGuides[iGuide].fFreq = pSpectralPeaks[iPeak].fFreq;
-			peakNotFound = 0;
-		}
-	}
-	return (1);
+    int iPeak = -1;
+    short peakNotFound = 1;
+
+    while (peakNotFound == 1 && *pFCurrentMax > 0)
+    {
+        /* \todo I don't think this ever returns -1, but check */
+        if ((iPeak = GetNextMax(pSpectralPeaks, pFCurrentMax)) < 0)
+            return (-1);
+
+        if (CheckForConflict (iPeak, pGuides, nGuides) < 0)
+        {
+            pGuides[iGuide].iPeakChosen = iPeak;
+            pGuides[iGuide].iStatus = GUIDE_BEG;
+            pGuides[iGuide].fFreq = pSpectralPeaks[iPeak].fFreq;
+            peakNotFound = 0;
+        }
+    }
+    return (1);
 }
 
 /*! \brief  function to advance the guides through the next frame
  *
  * the output is the frequency, magnitude, and phase tracks
  *
- * \param iFrame	 current frame number
+ * \param iFrame     current frame number
  * \param pAnalParams analysis parameters
  * \return error code \see SMS_ERRORS
  */
-int sms_peakContinuation (int iFrame, SMS_AnalParams *pAnalParams)
+int sms_peakContinuation(int iFrame, SMS_AnalParams *pAnalParams)
 {
-	int iGuide, iCurrentPeak = -1, iGoodPeak = -1;
-	sfloat fFund = pAnalParams->ppFrames[iFrame]->fFundamental,
-	       fFreqDev = fFund * pAnalParams->fFreqDeviation, fCurrentMax = 1000;
-  	static SMS_Guide *pGuides = NULL;
-  	static int nGuides = 0;
-
-  	if (pGuides == NULL || nGuides != pAnalParams->nGuides || pAnalParams->resetGuides == 1)
-	{
-		if ((pGuides = (SMS_Guide *) calloc(pAnalParams->nGuides, sizeof(SMS_Guide))) 
-			== NULL)
-			return (SMS_MALLOC);
-		if (pAnalParams->iFormat == SMS_FORMAT_H ||
-			pAnalParams->iFormat == SMS_FORMAT_HP)
-			for (iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
-				pGuides[iGuide].fFreq = pAnalParams->fDefaultFundamental
-				* (iGuide + 1);
-		nGuides = pAnalParams->nGuides;
-		pAnalParams->resetGuides = 0;
-	}
-
-	/* update guides with fundamental contribution */
-	if (fFund > 0 &&
-	    (pAnalParams->iFormat == SMS_FORMAT_H ||
-	     pAnalParams->iFormat == SMS_FORMAT_HP))
-		for(iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
-			pGuides[iGuide].fFreq = 
-				(1 - pAnalParams->fFundContToGuide) * pGuides[iGuide].fFreq +        
-				pAnalParams->fFundContToGuide * fFund * (iGuide + 1);
-  
-	if (pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
-	    pAnalParams->iDebugMode == SMS_DBG_ALL)
-		fprintf (stdout, 
-		         "Frame %d Peak Continuation: \n", 
-		         pAnalParams->ppFrames[iFrame]->iFrameNum);
-
-	/* continue all guides */
-	for (iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
-	{
-		sfloat fPreviousFreq = 
-			pAnalParams->ppFrames[iFrame-1]->deterministic.pFSinFreq[iGuide];
-
-		/* get the guide value by upgrading the previous guide */
-		if (fPreviousFreq > 0)
-			pGuides[iGuide].fFreq =
-				(1 - pAnalParams->fPeakContToGuide) * pGuides[iGuide].fFreq +
-				pAnalParams->fPeakContToGuide * fPreviousFreq;
-   
-		if (pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
-		    pAnalParams->iDebugMode == SMS_DBG_ALL)
-			fprintf (stdout, "Guide %d:  freq %f, mag %f\n", 
-			         iGuide, pGuides[iGuide].fFreq, pGuides[iGuide].fMag);
-      
-		if (pGuides[iGuide].fFreq <= 0.0 ||
-		    pGuides[iGuide].fFreq > pAnalParams->fHighestFreq)
-		{
-			pGuides[iGuide].iStatus = GUIDE_DEAD;
-			pGuides[iGuide].fFreq = 0;
-			continue;
-		}
-
-		pGuides[iGuide].iPeakChosen = -1;
-    
-		if (pAnalParams->iFormat == SMS_FORMAT_IH ||
-		    pAnalParams->iFormat == SMS_FORMAT_IHP)
-			fFreqDev = pGuides[iGuide].fFreq * pAnalParams->fFreqDeviation;
-
-		/* get the best peak for the guide */
-		iGoodPeak = 
-			GetBestPeak(pGuides, iGuide, pAnalParams->ppFrames[iFrame]->pSpectralPeaks, 
-			            pAnalParams, fFreqDev);
-
-//		printf("%f\t->\t%f (", fPreviousFreq, pAnalParams->ppFrames[iFrame]->pSpectralPeaks[iGoodPeak].fFreq);
-//		int i;
-//		for(i = 0; i < pAnalParams->ppFrames[iFrame]->nPeaks; i++)
-//		{
-//			printf("%f ", pAnalParams->ppFrames[iFrame]->pSpectralPeaks[i].fFreq);
-//		}
-//		printf(")\n");
-	}
-  
-	/* try to find good peaks for the GUIDE_DEAD guides */
-	if (pAnalParams->iFormat == SMS_FORMAT_IH ||
-	    pAnalParams->iFormat == SMS_FORMAT_IHP)
-		for(iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
-		{
-			if (pGuides[iGuide].iStatus != GUIDE_DEAD)
-				continue; 
-	
-			if (GetStartingPeak (iGuide, pGuides, pAnalParams->nGuides, 
-			                     pAnalParams->ppFrames[iFrame]->pSpectralPeaks,
-			                     &fCurrentMax) == -1)
-				break;
-		}
-
-	/* save all the continuation values,
-	 * assume output tracks are already clear */
-	for (iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
-	{
-		if (pGuides[iGuide].iStatus == GUIDE_DEAD)
-			continue; 
-
-		if (pAnalParams->iFormat == SMS_FORMAT_IH ||
-		    pAnalParams->iFormat == SMS_FORMAT_IHP)
-		{
-			if (pGuides[iGuide].iStatus > 0 &&
-			    pGuides[iGuide].iPeakChosen == -1)
-			{ 
-				if(pGuides[iGuide].iStatus++ > pAnalParams->iMaxSleepingTime)
-				{
-					pGuides[iGuide].iStatus = GUIDE_DEAD;
-					pGuides[iGuide].fFreq = 0;
-					pGuides[iGuide].fMag = 0;
-					pGuides[iGuide].iPeakChosen = -1;	  	  
-	 			}
-				else
-					pGuides[iGuide].iStatus++;
-				continue;
-			}
-      
-			if (pGuides[iGuide].iStatus == GUIDE_ACTIVE &&
-			    pGuides[iGuide].iPeakChosen == -1)
-			{
-				pGuides[iGuide].iStatus = 1;
-				continue;
-			}
-		}
-
-		/* if good continuation peak found, save it */
-		if ((iCurrentPeak = pGuides[iGuide].iPeakChosen) >= 0)
-		{
-			pAnalParams->ppFrames[iFrame]->deterministic.pFSinFreq[iGuide] = 
-				pAnalParams->ppFrames[iFrame]->pSpectralPeaks[iCurrentPeak].fFreq;
-			pAnalParams->ppFrames[iFrame]->deterministic.pFSinAmp[iGuide] =
-				pAnalParams->ppFrames[iFrame]->pSpectralPeaks[iCurrentPeak].fMag;
-			pAnalParams->ppFrames[iFrame]->deterministic.pFSinPha[iGuide] = 
-				pAnalParams->ppFrames[iFrame]->pSpectralPeaks[iCurrentPeak].fPhase;
-
-			pGuides[iGuide].iStatus = GUIDE_ACTIVE;
-			pGuides[iGuide].iPeakChosen = -1;
-		}
-	}
-        return(SMS_OK);
+    int iGuide, iCurrentPeak = -1, iGoodPeak = -1;
+    sfloat fFund = pAnalParams->ppFrames[iFrame]->fFundamental;
+    sfloat fFreqDev = fFund * pAnalParams->fFreqDeviation, fCurrentMax = 1000;
+    static SMS_Guide *pGuides = NULL;
+    static int nGuides = 0;
+
+    if(pGuides == NULL || nGuides != pAnalParams->nGuides || pAnalParams->resetGuides == 1)
+    {
+        if((pGuides = (SMS_Guide *) calloc(pAnalParams->nGuides, sizeof(SMS_Guide))) 
+           == NULL)
+            return SMS_MALLOC;
+        if (pAnalParams->iFormat == SMS_FORMAT_H ||
+            pAnalParams->iFormat == SMS_FORMAT_HP)
+            for (iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
+                pGuides[iGuide].fFreq = pAnalParams->fDefaultFundamental * (iGuide + 1);
+        nGuides = pAnalParams->nGuides;
+        pAnalParams->resetGuides = 0;
+    }
+
+    /* update guides with fundamental contribution */
+    if(fFund > 0 &&
+       (pAnalParams->iFormat == SMS_FORMAT_H ||
+        pAnalParams->iFormat == SMS_FORMAT_HP))
+        for(iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
+            pGuides[iGuide].fFreq = 
+                (1 - pAnalParams->fFundContToGuide) * pGuides[iGuide].fFreq +        
+                 pAnalParams->fFundContToGuide * fFund * (iGuide + 1);
+
+    if (pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
+            pAnalParams->iDebugMode == SMS_DBG_ALL)
+        fprintf(stdout, 
+                "Frame %d Peak Continuation: \n", 
+                pAnalParams->ppFrames[iFrame]->iFrameNum);
+
+    /* continue all guides */
+    for(iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
+    {
+        sfloat fPreviousFreq = 
+            pAnalParams->ppFrames[iFrame-1]->deterministic.pFSinFreq[iGuide];
+
+        /* get the guide value by upgrading the previous guide */
+        if(fPreviousFreq > 0)
+            pGuides[iGuide].fFreq =
+                (1 - pAnalParams->fPeakContToGuide) * pGuides[iGuide].fFreq +
+                pAnalParams->fPeakContToGuide * fPreviousFreq;
+
+        if(pAnalParams->iDebugMode == SMS_DBG_PEAK_CONT ||
+           pAnalParams->iDebugMode == SMS_DBG_ALL)
+            fprintf(stdout, "Guide %d:  freq %f, mag %f\n", 
+                    iGuide, pGuides[iGuide].fFreq, pGuides[iGuide].fMag);
+
+        if(pGuides[iGuide].fFreq <= 0.0 ||
+           pGuides[iGuide].fFreq > pAnalParams->fHighestFreq)
+        {
+            pGuides[iGuide].iStatus = GUIDE_DEAD;
+            pGuides[iGuide].fFreq = 0;
+            continue;
+        }
+
+        pGuides[iGuide].iPeakChosen = -1;
+
+        if(pAnalParams->iFormat == SMS_FORMAT_IH ||
+           pAnalParams->iFormat == SMS_FORMAT_IHP)
+            fFreqDev = pGuides[iGuide].fFreq * pAnalParams->fFreqDeviation;
+
+        /* get the best peak for the guide */
+       GetBestPeak(pGuides, iGuide, pAnalParams->ppFrames[iFrame]->pSpectralPeaks, 
+                   pAnalParams, fFreqDev);
+    }
+
+    /* try to find good peaks for the GUIDE_DEAD guides */
+    if(pAnalParams->iFormat == SMS_FORMAT_IH ||
+       pAnalParams->iFormat == SMS_FORMAT_IHP)
+    {
+        for(iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
+        {
+            if(pGuides[iGuide].iStatus != GUIDE_DEAD)
+                continue; 
+
+            /* TODO: make this function use pAnalParams->maxPeaks instead
+             * of SMS_MAX_NPEAKS */
+            if(GetStartingPeak (iGuide, pGuides, pAnalParams->nGuides, 
+               pAnalParams->ppFrames[iFrame]->pSpectralPeaks,
+               &fCurrentMax) == -1)
+                break;
+        }
+    }
+
+    /* save all the continuation values,
+     * assume output tracks are already clear */
+    for(iGuide = 0; iGuide < pAnalParams->nGuides; iGuide++)
+    {
+        if(pGuides[iGuide].iStatus == GUIDE_DEAD)
+            continue; 
+
+        if(pAnalParams->iFormat == SMS_FORMAT_IH ||
+           pAnalParams->iFormat == SMS_FORMAT_IHP)
+        {
+            if(pGuides[iGuide].iStatus > 0 &&
+               pGuides[iGuide].iPeakChosen == -1)
+            { 
+                if(pGuides[iGuide].iStatus++ > pAnalParams->iMaxSleepingTime)
+                {
+                    pGuides[iGuide].iStatus = GUIDE_DEAD;
+                    pGuides[iGuide].fFreq = 0;
+                    pGuides[iGuide].fMag = 0;
+                    pGuides[iGuide].iPeakChosen = -1;         
+                }
+                else
+                    pGuides[iGuide].iStatus++;
+                continue;
+            }
+
+            if(pGuides[iGuide].iStatus == GUIDE_ACTIVE &&
+               pGuides[iGuide].iPeakChosen == -1)
+            {
+                pGuides[iGuide].iStatus = 1;
+                continue;
+            }
+        }
+
+        /* if good continuation peak found, save it */
+        if((iCurrentPeak = pGuides[iGuide].iPeakChosen) >= 0)
+        {
+            pAnalParams->ppFrames[iFrame]->deterministic.pFSinFreq[iGuide] = 
+                pAnalParams->ppFrames[iFrame]->pSpectralPeaks[iCurrentPeak].fFreq;
+            pAnalParams->ppFrames[iFrame]->deterministic.pFSinAmp[iGuide] =
+                pAnalParams->ppFrames[iFrame]->pSpectralPeaks[iCurrentPeak].fMag;
+            pAnalParams->ppFrames[iFrame]->deterministic.pFSinPha[iGuide] = 
+                pAnalParams->ppFrames[iFrame]->pSpectralPeaks[iCurrentPeak].fPhase;
+
+            pGuides[iGuide].iStatus = GUIDE_ACTIVE;
+            pGuides[iGuide].iPeakChosen = -1;
+        }
+    }
+    return SMS_OK;
 }
diff --git a/tests/sms.py b/tests/sms.py
index 1159620..36cc875 100644
--- a/tests/sms.py
+++ b/tests/sms.py
@@ -22,11 +22,11 @@ from scipy.io.wavfile import read
 import unittest
 
 class TestSimplSMS(unittest.TestCase):
-    FLOAT_PRECISION = 5 # number of decimal places to check for accuracy
+    FLOAT_PRECISION = 3 # number of decimal places to check for accuracy
     input_file = 'audio/flute.wav'
     frame_size = 2048
     hop_size = 512
-    num_frames = 20
+    num_frames = 9 
     num_samples = frame_size + ((num_frames - 1) * hop_size)
     max_peaks = 10
     max_partials = 3
@@ -50,12 +50,14 @@ class TestSimplSMS(unittest.TestCase):
         analysis_params.iWindowType = pysms.SMS_WIN_HAMMING
         analysis_params.fDefaultFundamental = 100
         analysis_params.fHighestFreq = 20000
+        analysis_params.iFormat = pysms.SMS_FORMAT_HP
         analysis_params.nTracks = max_peaks
-        analysis_params.maxPeaks = max_peaks
-        analysis_params.iMaxDelayFrames = num_frames + 1
+        #analysis_params.nGuides = max_peaks
+        analysis_params.iMaxDelayFrames = 4
         analysis_params.analDelay = 0
         analysis_params.minGoodFrames = 1
         analysis_params.iCleanTracks = 0
+        analysis_params.iStochasticType = pysms.SMS_STOC_NONE
         pysms.sms_initAnalysis(analysis_params, snd_header)
         sms_header.nStochasticCoeff = 128
         pysms.sms_fillHeader(sms_header, analysis_params, "pysms")
@@ -70,18 +72,21 @@ class TestSimplSMS(unittest.TestCase):
         analysis_params.iWindowType = simplsms.SMS_WIN_HAMMING
         analysis_params.fDefaultFundamental = 100
         analysis_params.fHighestFreq = 20000
-        analysis_params.iMaxDelayFrames = num_frames + 1
+        analysis_params.iFormat = simplsms.SMS_FORMAT_HP
+        analysis_params.nTracks = max_peaks
+        analysis_params.maxPeaks = max_peaks
+        #analysis_params.nGuides = max_peaks
+        analysis_params.resetGuides = 1
+        analysis_params.iMaxDelayFrames = 4
         analysis_params.analDelay = 0
         analysis_params.minGoodFrames = 1
         analysis_params.iCleanTracks = 0
-        analysis_params.iFormat = pysms.SMS_FORMAT_HP
-        analysis_params.nTracks = max_peaks
-        analysis_params.maxPeaks = max_peaks
+        analysis_params.iStochasticType = simplsms.SMS_STOC_NONE
         simplsms.sms_initAnalysis(analysis_params)
         sms_header.nStochasticCoeff = 128
         simplsms.sms_fillHeader(sms_header, analysis_params, "simplsms")
         return analysis_params, sms_header
-       
+
     def test_size_next_read(self):
         """test_size_next_read
         Make sure pysms PeakDetection is calculating 
@@ -90,11 +95,14 @@ class TestSimplSMS(unittest.TestCase):
         analysis_params, sms_header, snd_header = self.pysms_params(sampling_rate,
                                                                     self.max_peaks,
                                                                     self.num_frames)
+        analysis_params.iMaxDelayFrames = self.num_frames + 1
+        pysms.sms_initAnalysis(analysis_params, snd_header)
+
         sample_offset = 0
         pysms_size_new_data = 0
         current_frame = 0
         sms_next_read_sizes = []
-        
+
         while current_frame < self.num_frames:
             sms_next_read_sizes.append(analysis_params.sizeNextRead)
             sample_offset += pysms_size_new_data
@@ -117,19 +125,19 @@ class TestSimplSMS(unittest.TestCase):
         pysms.sms_freeAnalysis(analysis_params)
         pysms.sms_closeSF()
         pysms.sms_free()
-        
+
         pd = simpl.SMSPeakDetection()
         pd.hop_size = self.hop_size
         current_frame = 0
         sample_offset = 0
-        
+
         while current_frame < self.num_frames:
             pd.frame_size = pd.get_next_frame_size()
             self.assertEquals(sms_next_read_sizes[current_frame], pd.frame_size)
             pd.find_peaks_in_frame(audio[sample_offset:sample_offset + pd.frame_size])
             sample_offset += pd.frame_size
             current_frame += 1
-            
+
     def test_peak_detection(self): 
         """test_peak_detection
         Compare simplsms Peaks with SMS peaks. Exact peak
@@ -145,11 +153,14 @@ class TestSimplSMS(unittest.TestCase):
         analysis_params, sms_header = self.simplsms_params(sampling_rate,
                                                            self.max_peaks,
                                                            self.num_frames)
+        analysis_params.iMaxDelayFrames = self.num_frames + 1
+        simplsms.sms_initAnalysis(analysis_params)
+
         sample_offset = 0
         size_new_data = 0
         current_frame = 0
         sms_peaks = []
- 
+
         while current_frame < self.num_frames:
             sample_offset += size_new_data
             size_new_data = analysis_params.sizeNextRead
@@ -194,11 +205,11 @@ class TestSimplSMS(unittest.TestCase):
         while current_frame < self.num_frames:
             pd.frame_size = pd.get_next_frame_size()
             simpl_peaks.append(
-                pd.find_peaks_in_frame(audio[sample_offset:sample_offset + pd.frame_size]))
+                    pd.find_peaks_in_frame(audio[sample_offset:sample_offset + pd.frame_size]))
             sample_offset += pd.frame_size
             #print current_frame, self.num_samples, pd.frame_size, sample_offset
             current_frame += 1
-            
+
         # make sure we have the same number of frames
         self.assertEquals(len(sms_peaks), len(simpl_peaks))
 
@@ -214,120 +225,307 @@ class TestSimplSMS(unittest.TestCase):
                 sms_peak = sms_frame[peak_number]
                 simpl_peak = simpl_frame[peak_number]
                 self.assertAlmostEquals(sms_peak.amplitude, simpl_peak.amplitude,
-                                       places=self.FLOAT_PRECISION)
+                                        places=self.FLOAT_PRECISION)
                 self.assertAlmostEquals(sms_peak.frequency, simpl_peak.frequency,
-                                       places=self.FLOAT_PRECISION)
+                                        places=self.FLOAT_PRECISION)
                 self.assertAlmostEquals(sms_peak.phase, simpl_peak.phase,
-                                       places=self.FLOAT_PRECISION)  
+                                        places=self.FLOAT_PRECISION)  
 
-    #def test_partial_tracking(self):
-    #    """test_partial_tracking
-    #    Compare pysms Partials with SMS partials.""" 
-    #    pysms.sms_init()
-    #    sms_header = pysms.SMS_Header()
-    #    snd_header = pysms.SMS_SndHeader()
-    #    # Try to open the input file to fill snd_header
-    #    if(pysms.sms_openSF(input_file, snd_header)):
-    #        raise NameError("error opening sound file: " + pysms.sms_errorString())
-    #    analysis_params = pysms.SMS_AnalParams()
-    #    analysis_params.iSamplingRate = sampling_rate
-    #    analysis_params.iFrameRate = sampling_rate / hop_size
-    #    sms_header.nStochasticCoeff = 128
-    #    analysis_params.fDefaultFundamental = 100
-    #    analysis_params.fHighestFreq = 20000
-    #    analysis_params.iMaxDelayFrames = 3
-    #    analysis_params.analDelay = 0
-    #    analysis_params.minGoodFrames = 1
-    #    analysis_params.iFormat = pysms.SMS_FORMAT_HP
-    #    analysis_params.nTracks = max_partials
-    #    analysis_params.nGuides = max_partials
-    #    analysis_params.iWindowType = pysms.SMS_WIN_HAMMING
-    #    pysms.sms_initAnalysis(analysis_params, snd_header)
-    #    analysis_params.nFrames = num_samples / hop_size
-    #    analysis_params.iSizeSound = num_samples
-    #    analysis_params.peakParams.iMaxPeaks = max_peaks
-    #    analysis_params.iStochasticType = pysms.SMS_STOC_NONE
-    #    pysms.sms_fillHeader(sms_header, analysis_params, "pysms")
+    def test_sms_analyze(self):
+        """test_sms_analyze
+        Make sure that the simplsms.sms_analyze function does the same thing
+        as the sms_analyze function from libsms."""
+        audio, sampling_rate = self.get_audio()
+        analysis_params, sms_header, snd_header = self.pysms_params(sampling_rate,
+                                                                    self.max_peaks,
+                                                                    self.num_frames)
+        analysis_params.nFrames = self.num_frames
+        analysis_params.iSizeSound = self.num_samples
+        analysis_params.peakParams.iMaxPeaks = self.max_peaks
 
-    #    sample_offset = 0
-    #    size_new_data = 0
-    #    current_frame = 0
-    #    sms_partials = []
-    #    live_partials = [None for i in range(max_partials)]
-    #    do_analysis = True
+        sample_offset = 0
+        size_new_data = 0
+        current_frame = 0
+        sms_partials = []
+        live_partials = [None for i in range(self.max_peaks)]
+        do_analysis = True
 
-    #    while do_analysis and (current_frame < num_frames):
-    #        sample_offset += size_new_data
-    #        if((sample_offset + analysis_params.sizeNextRead) < num_samples):
-    #            size_new_data = analysis_params.sizeNextRead
-    #        else:
-    #            size_new_data = num_samples - sample_offset
-    #        frame = audio[sample_offset:sample_offset + size_new_data]
-    #        analysis_data = pysms.SMS_Data()
-    #        pysms.sms_allocFrameH(sms_header, analysis_data)
-    #        status = pysms.sms_analyze(frame, analysis_data, analysis_params)  
-            
-    #        if status == 1:
-    #            num_partials = analysis_data.nTracks
-    #            pysms_freqs = simpl.zeros(num_partials)
-    #            pysms_amps = simpl.zeros(num_partials)
-    #            pysms_phases = simpl.zeros(num_partials)
-    #            analysis_data.getSinFreq(pysms_freqs)
-    #            analysis_data.getSinAmp(pysms_amps)
-    #            analysis_data.getSinPhase(pysms_phases)
-    #            # make partial objects
-    #            for i in range(num_partials):
-    #                # for each partial, if the mag is > 0, this partial is alive
-    #                if pysms_amps[i] > 0:
-    #                    # create a peak object
-    #                    p = simpl.Peak()
-    #                    p.amplitude = pysms_amps[i]
-    #                    p.frequency = pysms_freqs[i]
-    #                    p.phase = pysms_phases[i]
-    #                    # add this peak to the appropriate partial
-    #                    if not live_partials[i]:
-    #                        live_partials[i] = simpl.Partial()
-    #                        live_partials[i].starting_frame = current_frame
-    #                        sms_partials.append(live_partials[i])
-    #                    live_partials[i].add_peak(p)
-    #                # if the mag is 0 and this partial was alive, kill it
-    #                else:
-    #                    if live_partials[i]:
-    #                        live_partials[i] = None
-    #        elif status == -1:
-    #            do_analysis = False
-    #        current_frame += 1
+        while do_analysis and (current_frame < self.num_frames):
+            sample_offset += size_new_data
+            size_new_data = analysis_params.sizeNextRead
+            # convert frame to floats for libsms
+            frame = audio[sample_offset:sample_offset + size_new_data]
+            frame = np.array(frame, dtype=np.float32)
+            analysis_data = pysms.SMS_Data()
+            pysms.sms_allocFrameH(sms_header, analysis_data)
+            status = pysms.sms_analyze(frame, analysis_data, analysis_params)  
+
+            if status == 1:
+                num_partials = analysis_data.nTracks
+                sms_freqs = np.zeros(num_partials, dtype=np.float32)
+                sms_amps = np.zeros(num_partials, dtype=np.float32)
+                sms_phases = np.zeros(num_partials, dtype=np.float32)
+                analysis_data.getSinFreq(sms_freqs)
+                analysis_data.getSinAmp(sms_amps)
+                analysis_data.getSinPhase(sms_phases)
+                # make partial objects
+                for i in range(num_partials):
+                    # for each partial, if the mag is > 0, this partial is alive
+                    if sms_amps[i] > 0:
+                        # create a peak object
+                        p = simpl.Peak()
+                        p.amplitude = sms_amps[i]
+                        p.frequency = sms_freqs[i]
+                        p.phase = sms_phases[i]
+                        # add this peak to the appropriate partial
+                        if not live_partials[i]:
+                            live_partials[i] = simpl.Partial()
+                            live_partials[i].starting_frame = current_frame
+                            sms_partials.append(live_partials[i])
+                        live_partials[i].add_peak(p)
+                    # if the mag is 0 and this partial was alive, kill it
+                else:
+                    if live_partials[i]:
+                        live_partials[i] = None
+            elif status == -1:
+                do_analysis = False
+            current_frame += 1
+
+        pysms.sms_freeAnalysis(analysis_params)
+        pysms.sms_closeSF()
+        pysms.sms_free()
+
+        return
+        simpl_analysis_params, simpl_sms_header = self.simplsms_params(sampling_rate,
+                                                                       self.max_peaks,
+                                                                       self.num_frames)
+        simpl_analysis_params.nFrames = self.num_frames
+        simpl_analysis_params.iSizeSound = self.num_samples
+
+        sample_offset = 0
+        size_new_data = 0
+        current_frame = 0
+        simplsms_partials = []
+        live_partials = [None for i in range(self.max_peaks)]
+        do_analysis = True
+
+        while do_analysis and (current_frame < self.num_frames):
+            sample_offset += size_new_data
+            size_new_data = simpl_analysis_params.sizeNextRead
+            frame = audio[sample_offset:sample_offset + size_new_data]
+            analysis_data = simplsms.SMS_Data()
+            simplsms.sms_allocFrameH(simpl_sms_header, analysis_data)
+            status = simplsms.sms_analyze(frame, analysis_data, simpl_analysis_params)  
+
+            if status == 1:
+                num_partials = analysis_data.nTracks
+                freqs = simpl.zeros(num_partials)
+                amps = simpl.zeros(num_partials)
+                phases = simpl.zeros(num_partials)
+                analysis_data.getSinFreq(freqs)
+                analysis_data.getSinAmp(amps)
+                analysis_data.getSinPhase(phases)
+                # make partial objects
+                for i in range(num_partials):
+                    # for each partial, if the mag is > 0, this partial is alive
+                    if amps[i] > 0:
+                        # create a peak object
+                        p = simpl.Peak()
+                        p.amplitude = amps[i]
+                        p.frequency = freqs[i]
+                        p.phase = phases[i]
+                        # add this peak to the appropriate partial
+                        if not live_partials[i]:
+                            live_partials[i] = simpl.Partial()
+                            live_partials[i].starting_frame = current_frame
+                            simplsms_partials.append(live_partials[i])
+                        live_partials[i].add_peak(p)
+                    # if the mag is 0 and this partial was alive, kill it
+                else:
+                    if live_partials[i]:
+                        live_partials[i] = None
+            elif status == -1:
+                do_analysis = False
+            current_frame += 1
+
+        simplsms.sms_freeAnalysis(simpl_analysis_params)
+        simplsms.sms_free()
+
+        return
+
+        # make sure both have the same number of partials
+        self.assertEquals(len(sms_partials), len(partials))
+
+        # make sure each partial is the same
+        for i in range(len(sms_partials)):
+            self.assertEquals(sms_partials[i].get_length(), partials[i].get_length())
+            for peak_number in range(sms_partials[i].get_length()):
+                self.assertAlmostEquals(sms_partials[i].peaks[peak_number].amplitude,
+                                        partials[i].peaks[peak_number].amplitude,
+                                        places = self.FLOAT_PRECISION)
+                self.assertAlmostEquals(sms_partials[i].peaks[peak_number].frequency,
+                                        partials[i].peaks[peak_number].frequency,
+                                        places = self.FLOAT_PRECISION)
+                self.assertAlmostEquals(sms_partials[i].peaks[peak_number].phase,
+                                        partials[i].peaks[peak_number].phase,
+                                        places = self.FLOAT_PRECISION)
+
+    def test_partial_tracking(self):
+        """test_partial_tracking
+        Compare pysms Partials with SMS partials.""" 
+        audio, sampling_rate = self.get_audio()
+        analysis_params, sms_header, snd_header = self.pysms_params(sampling_rate,
+                                                                    self.max_peaks,
+                                                                    self.num_frames)
+        analysis_params.nFrames = self.num_frames
+        analysis_params.iSizeSound = self.num_samples
+        analysis_params.peakParams.iMaxPeaks = self.max_peaks
+
+        sample_offset = 0
+        size_new_data = 0
+        current_frame = 0
+        sms_partials = []
+        live_partials = [None for i in range(self.max_peaks)]
+        do_analysis = True
+
+        while do_analysis and (current_frame < self.num_frames):
+            sample_offset += size_new_data
+            size_new_data = analysis_params.sizeNextRead
+            # convert frame to floats for libsms
+            frame = audio[sample_offset:sample_offset + size_new_data]
+            frame = np.array(frame, dtype=np.float32)
+            analysis_data = pysms.SMS_Data()
+            pysms.sms_allocFrameH(sms_header, analysis_data)
+            status = pysms.sms_analyze(frame, analysis_data, analysis_params)  
+
+            if status == 1:
+                num_partials = analysis_data.nTracks
+                sms_freqs = np.zeros(num_partials, dtype=np.float32)
+                sms_amps = np.zeros(num_partials, dtype=np.float32)
+                sms_phases = np.zeros(num_partials, dtype=np.float32)
+                analysis_data.getSinFreq(sms_freqs)
+                analysis_data.getSinAmp(sms_amps)
+                analysis_data.getSinPhase(sms_phases)
+                # make partial objects
+                for i in range(num_partials):
+                    # for each partial, if the mag is > 0, this partial is alive
+                    if sms_amps[i] > 0:
+                        # create a peak object
+                        p = simpl.Peak()
+                        p.amplitude = sms_amps[i]
+                        p.frequency = sms_freqs[i]
+                        p.phase = sms_phases[i]
+                        # add this peak to the appropriate partial
+                        if not live_partials[i]:
+                            live_partials[i] = simpl.Partial()
+                            live_partials[i].starting_frame = current_frame
+                            sms_partials.append(live_partials[i])
+                        live_partials[i].add_peak(p)
+                    # if the mag is 0 and this partial was alive, kill it
+                else:
+                    if live_partials[i]:
+                        live_partials[i] = None
+            elif status == -1:
+                do_analysis = False
+            current_frame += 1
+
+        pysms.sms_freeAnalysis(analysis_params)
+        pysms.sms_closeSF()
+        pysms.sms_free()
+
+        import debug
         
-    #    pysms.sms_freeAnalysis(analysis_params)
-    #    pysms.sms_closeSF()
-    #    pysms.sms_free()
+        debug.print_partials(sms_partials)
+        print
 
-    #    pd = simpl.SMSPeakDetection()
-    #    pd.max_peaks = max_peaks
-    #    pd.hop_size = hop_size 
-    #    peaks = pd.find_peaks(audio)
-    #    pt = simpl.SMSPartialTracking()
-    #    pt.max_partials = max_partials
-    #    partials = pt.find_partials(peaks[0:num_frames])
-               
-    #    # make sure both have the same number of partials
-    #    self.assertEquals(len(sms_partials), len(partials))
-
-    #    # make sure each partial is the same
-    #    for i in range(len(sms_partials)):
-    #        self.assertEquals(sms_partials[i].get_length(), partials[i].get_length())
-    #        for peak_number in range(sms_partials[i].get_length()):
-    #            self.assertAlmostEquals(sms_partials[i].peaks[peak_number].amplitude,
-    #                                    partials[i].peaks[peak_number].amplitude,
-    #                                    places = FLOAT_PRECISION)
-    #            self.assertAlmostEquals(sms_partials[i].peaks[peak_number].frequency,
-    #                                    partials[i].peaks[peak_number].frequency,
-    #                                    places = FLOAT_PRECISION)
-    #            self.assertAlmostEquals(sms_partials[i].peaks[peak_number].phase,
-    #                                    partials[i].peaks[peak_number].phase,
-    #                                    places = FLOAT_PRECISION)
-                
-    #def test_interpolate_frames(self):
+        analysis_params, sms_header = self.simplsms_params(sampling_rate,
+                                                           self.max_peaks,
+                                                           self.num_frames)
+        analysis_params.nFrames = self.num_frames
+        analysis_params.iSizeSound = self.num_samples
+
+        sample_offset = 0
+        size_new_data = 0
+        current_frame = 0
+        sms_partials = []
+        live_partials = [None for i in range(self.max_peaks)]
+        do_analysis = True
+
+        while do_analysis and (current_frame < self.num_frames):
+            sample_offset += size_new_data
+            size_new_data = analysis_params.sizeNextRead
+            frame = audio[sample_offset:sample_offset + size_new_data]
+            analysis_data = simplsms.SMS_Data()
+            simplsms.sms_allocFrameH(sms_header, analysis_data)
+            status = simplsms.sms_analyze(frame, analysis_data, analysis_params)  
+
+            if status == 1:
+                num_partials = analysis_data.nTracks
+                sms_freqs = np.zeros(num_partials, dtype=np.float64)
+                sms_amps = np.zeros(num_partials, dtype=np.float64)
+                sms_phases = np.zeros(num_partials, dtype=np.float64)
+                analysis_data.getSinFreq(sms_freqs)
+                analysis_data.getSinAmp(sms_amps)
+                analysis_data.getSinPhase(sms_phases)
+                # make partial objects
+                for i in range(num_partials):
+                    # for each partial, if the mag is > 0, this partial is alive
+                    if sms_amps[i] > 0:
+                        # create a peak object
+                        p = simpl.Peak()
+                        p.amplitude = sms_amps[i]
+                        p.frequency = sms_freqs[i]
+                        p.phase = sms_phases[i]
+                        # add this peak to the appropriate partial
+                        if not live_partials[i]:
+                            live_partials[i] = simpl.Partial()
+                            live_partials[i].starting_frame = current_frame
+                            sms_partials.append(live_partials[i])
+                        live_partials[i].add_peak(p)
+                    # if the mag is 0 and this partial was alive, kill it
+                else:
+                    if live_partials[i]:
+                        live_partials[i] = None
+            elif status == -1:
+                do_analysis = False
+            current_frame += 1
+
+        simplsms.sms_freeAnalysis(analysis_params)
+        simplsms.sms_free()
+
+        debug.print_partials(sms_partials)
+        return
+
+        #pd = simpl.SMSPeakDetection()
+        #pd.max_peaks = self.max_peaks
+        #pd.hop_size = self.hop_size 
+        #peaks = pd.find_peaks(audio)
+        #pt = simpl.SMSPartialTracking()
+        #pt.max_partials = self.max_peaks
+        #partials = pt.find_partials(peaks[0:self.num_frames])
+
+        #from pylab import show
+        #simpl.plot.plot_partials(sms_partials)
+        #show()
+
+        # make sure both have the same number of partials
+        self.assertEquals(len(sms_partials), len(partials))
+
+        # make sure each partial is the same
+        for i in range(len(sms_partials)):
+            self.assertEquals(sms_partials[i].get_length(), partials[i].get_length())
+            for peak_number in range(sms_partials[i].get_length()):
+                self.assertAlmostEquals(sms_partials[i].peaks[peak_number].amplitude,
+                                        partials[i].peaks[peak_number].amplitude,
+                                        places = self.FLOAT_PRECISION)
+                self.assertAlmostEquals(sms_partials[i].peaks[peak_number].frequency,
+                                        partials[i].peaks[peak_number].frequency,
+                                        places = self.FLOAT_PRECISION)
+                self.assertAlmostEquals(sms_partials[i].peaks[peak_number].phase,
+                                        partials[i].peaks[peak_number].phase,
+                                        places = self.FLOAT_PRECISION)
+
+               #def test_interpolate_frames(self):
     #    """test_interpolate_frames
     #    Make sure that pysms.sms_interpolateFrames returns the expected values
     #    with interpolation factors of 0 and 1."""
@@ -424,10 +622,10 @@ class TestSimplSMS(unittest.TestCase):
     #        elif status == -1:
     #            raise Exception("AnalysisStoppedEarly")
     #        current_frame += 1
-        
+
     #    pysms.sms_freeAnalysis(analysis_params)
     #    pysms.sms_closeSF()
-        
+
     #def test_harmonic_synthesis(self):
     #    """test_harmonic_synthesis
     #    Compare pysms synthesised harmonic component with SMS synthesised 
@@ -479,10 +677,10 @@ class TestSimplSMS(unittest.TestCase):
     #        if status == -1:
     #            do_analysis = False
     #        current_frame += 1
-        
+
     #    pysms.sms_freeAnalysis(analysis_params)
     #    pysms.sms_closeSF()
-        
+
     #    interp_frame = pysms.SMS_Data() 
     #    synth_params = pysms.SMS_SynthParams() 
     #    synth_params.iSynthesisType = pysms.SMS_STYPE_DET
@@ -504,10 +702,10 @@ class TestSimplSMS(unittest.TestCase):
     #        pysms_audio = np.hstack((pysms_audio, synth_samples))
     #        num_synth_samples += synth_params.sizeHop
     #        current_frame += 1
-              
+
     #    pysms.sms_freeSynth(synth_params)
     #    pysms.sms_free()     
-        
+
     #    pd = simpl.SMSPeakDetection()
     #    pd.max_peaks = max_peaks
     #    pd.hop_size = hop_size 
@@ -526,7 +724,7 @@ class TestSimplSMS(unittest.TestCase):
     #    for i in range(simpl_audio.size):
     #        self.assertAlmostEquals(pysms_audio[i], simpl_audio[i],
     #                                places = FLOAT_PRECISION)
-            
+
     #def test_residual_synthesis(self):
     #    """test_residual_synthesis
     #    Compare pysms residual signal with SMS residual""" 
@@ -577,11 +775,11 @@ class TestSimplSMS(unittest.TestCase):
     #        if status == -1:
     #            do_analysis = False
     #        current_frame += 1
-        
+
     #    pysms.sms_freeAnalysis(analysis_params)
     #    pysms.sms_closeSF()
     #    pysms.sms_free()
-        
+
     #    pd = simpl.SMSPeakDetection()
     #    pd.max_peaks = max_peaks
     #    pd.hop_size = hop_size
@@ -617,8 +815,9 @@ class TestSimplSMS(unittest.TestCase):
 if __name__ == '__main__':
     suite = unittest.TestSuite()
     suite.addTest(TestSimplSMS('test_size_next_read'))
-    suite.addTest(TestSimplSMS('test_peak_detection'))
-    #suite.addTest(TestSimplSMS('test_partial_tracking'))
+    #suite.addTest(TestSimplSMS('test_peak_detection'))
+    #suite.addTest(TestSimplSMS('test_sms_analyze'))
+    suite.addTest(TestSimplSMS('test_partial_tracking'))
     #suite.addTest(TestSimplSMS('test_interpolate_frames'))
     #suite.addTest(TestSimplSMS('test_harmonic_synthesis'))
     #suite.addTest(TestSimplSMS('test_residual_synthesis'))