summaryrefslogtreecommitdiff
path: root/sms/peakContinuation.c
diff options
context:
space:
mode:
Diffstat (limited to 'sms/peakContinuation.c')
-rw-r--r--sms/peakContinuation.c513
1 files changed, 513 insertions, 0 deletions
diff --git a/sms/peakContinuation.c b/sms/peakContinuation.c
new file mode 100644
index 0000000..7af4388
--- /dev/null
+++ b/sms/peakContinuation.c
@@ -0,0 +1,513 @@
+/*
+ * Copyright (c) 2008 MUSIC TECHNOLOGY GROUP (MTG)
+ * UNIVERSITAT POMPEU FABRA
+ *
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+ *
+ */
+/*! \file peakContinuation.c
+ * \brief peak continuation algorithm and functions
+ */
+
+#include "sms.h"
+
+/*! diferent status of guide */
+#define GUIDE_BEG -2
+#define GUIDE_DEAD -1
+#define GUIDE_ACTIVE 0
+
+#define MAX_CONT_CANDIDATES 5 /*!< maximum number of peak continuation
+ 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
+ * \return peak number or -1 if nothing is good
+ */
+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);
+}
+
+/*! \brief choose the best candidate out of all
+ *
+ * \param pCandidate pointer to all the continuation candidates
+ * \param nCandidates number of candidates
+ * \param fFreqDev maximum frequency deviation allowed
+ * \return the peak number of the best candidate
+ */
+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);
+}
+
+/*! \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
+ * \return number of guide that chose the peak, or -1 if none
+ */
+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;
+}
+
+/*! \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
+ * \return number of guide
+ */
+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);
+}
+
+/*! \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 fFreqDev frequency deviation allowed
+ * \return the peak number
+ */
+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);
+}
+
+/*! \brief function to get the next maximum (magnitude) peak
+ * \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);
+}
+
+/*! \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
+ * \return \todo should this return something?
+ */
+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);
+}
+
+/*! \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 pAnalParams analysis parameters
+ * \return error code \see SMS_ERRORS
+ */
+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);
+}