stft_highres.c

LastWave

				   STFT stftReal,STFT stftPhase,
				   STFT subStftReal,STFT subStftPhase)
{
    /* Checking arguments */
    CheckStftHighRes(stftHighRes);
    CheckStftReal(stftReal);
    CheckStftPhase(stftPhase);
    CheckStftReal(subStftReal);
    CheckStftPhase(subStftPhase);
    if(stftReal->flagUpToDate == NO)
	Errorf("CheckStftHighResEnergy : stftReal is out of date");
    if(stftPhase->flagUpToDate == NO)
	Errorf("CheckStftHighResEnergy : stftPhase is out of date");
    if(subStftReal->flagUpToDate == NO)
	Errorf("CheckStftHighResEnergy : subStftReal is out of date");
    if(subStftPhase->flagUpToDate == NO)
	Errorf("CheckStftHighResEnergy : subStftPhase is out of date");

    /* Checking coherence of time-freq grids */
    CheckTFContentCompat(stftReal,stftHighRes);
    CheckTFContentCompat(stftPhase,stftHighRes);
    CheckTFContentCompat(subStftReal,stftHighRes);
    CheckTFContentCompat(subStftPhase,stftHighRes);
    if(stftReal->windowSize        != stftHighRes->windowSize ||
       stftPhase->windowSize        != stftHighRes->windowSize ||
       subStftReal->windowSize     != subStftPhase->windowSize)
	Errorf("CheckStftHighResEnergy : incompatible windowSize");
    if(stftReal->borderType     != stftHighRes->borderType  ||  
       stftPhase->borderType     != stftHighRes->borderType  ||  
       subStftReal->borderType  != stftHighRes->borderType  ||  
       subStftPhase->borderType  != stftHighRes->borderType  ||  
       stftReal->windowShape    != stftHighRes->windowShape ||
       stftPhase->windowShape    != stftHighRes->windowShape ||
       subStftReal->windowShape != stftHighRes->windowShape ||
       subStftPhase->windowShape != stftHighRes->windowShape)
	Errorf("CheckStftHighResEnergy : incompatible borderType or windowShape");

    /* Later to remove */
    if(subStftReal->windowSize > stftReal->windowSize)
      Errorf("CheckStftHighResEnergy : subStft windowSize %d >= %d !",subStftReal->windowSize,stftReal->windowSize);

    /* Remark : (sub)RateEnergy should ALWAYS divide rateHighRes */
    if(stftHighRes->fRate % stftReal->fRate != 0)
      Errorf("CheckStftHighResEnergy : incompatible fRates %d %d",stftHighRes->fRate,stftReal->fRate);
    if(stftHighRes->tRate % stftReal->tRate != 0)
      Errorf("CheckStftHighResEnergy : incompatible tRates %d %d",stftHighRes->tRate,stftReal->tRate);
    if(stftHighRes->fRate % subStftReal->fRate != 0)
      Errorf("CheckStftHighResEnergy : incompatible sub fRates %d %d",stftHighRes->fRate,subStftReal->fRate);
    if(stftHighRes->tRate % subStftReal->tRate != 0)
      Errorf("CheckStftHighResEnergy : incompatible sub tRates %d %d",stftHighRes->tRate,subStftReal->tRate);

    if(stftHighRes->fRate % stftPhase->fRate != 0)
      Errorf("CheckStftHighResEnergy : incompatible fRates %d %d",stftHighRes->fRate,stftPhase->fRate);
    if(stftHighRes->tRate % stftPhase->tRate != 0)
      Errorf("CheckStftHighResEnergy : incompatible tRates %d %d",stftHighRes->tRate,stftPhase->tRate);
    if(stftHighRes->fRate % subStftPhase->fRate != 0)
      Errorf("CheckStftHighResEnergy : incompatible sub fRates %d %d",stftHighRes->fRate,subStftPhase->fRate);
    if(stftHighRes->tRate % subStftPhase->tRate != 0)
      Errorf("CheckStftHighResEnergy : incompatible sub tRates %d %d",stftHighRes->tRate,subStftPhase->tRate);
}



/***************************************************/
/*
 * Very helpful function to set a stft REAL HIGHRES part
 * at a given time 'timeId', for all frequencies.
 */
/***************************************************/
static void StftHighResSetFreq(STFT stftHighRes,
			       STFT stftReal,STFT stftPhase,
			       STFT subStftReal,STFT subStftPhase,
			       int timeId)
{
  static ATOM atom = NULL;
  int freqId;
  LWFLOAT *coeff2s,*phases,*coeff2sHighRes;
  LWFLOAT phase;
  /* Checkings */
  CheckStftHighResEnergy(stftHighRes,stftReal,stftPhase,subStftReal,subStftPhase);
  if(!INRANGE(0,timeId,stftHighRes->signalSize-1))
    Errorf("StftHighResSetFreq : bad timeId %d not in [0 %d]",timeId,stftHighRes->signalSize-1);

  /* Initializing (just once) */
  if(atom == NULL) {
    atom = NewAtom();
  } else {
    ClearAtom(atom);
  }
  CopyFieldsTFContent(stftHighRes,atom);
  atom->windowShape	= stftHighRes->windowShape;
  atom->windowSize = stftHighRes->windowSize;


  atom->timeId 	= timeId;
  atom->freqId  = 0;

  /* Getting the pointers to the data */
  GetStftData(stftReal,timeId,&coeff2s,NULL);
  GetStftData(stftHighRes,timeId,&coeff2sHighRes,NULL);
  GetStftData(stftPhase,timeId,&phases,NULL);

  /* Loop on all frequencies */
  for(freqId = 0;
      freqId < stftHighRes->fRate*stftHighRes->nSubBands; 
      freqId += stftHighRes->fRate) {
    // Setting the complex part of the atom
    atom->freqId  	 = freqId;
    // These fields are not used here
    atom->coeffR = 0.0;
    atom->coeffI = 0.0;
    // Convert to a real atom
    // Getting the energy and phase from the stftReal data
    if(!ComputeWindowGG(atom->windowShape,atom->windowSize,atom->freqId,&(atom->realGG),&(atom->imagGG)))
	Errorf("StftHighResSetFreq : error with ComputeWindowGG");

    atom->coeff2     = coeff2s[freqId/stftReal->fRate];
    phase            =   phases[freqId/stftReal->fRate];
    // DEBUG : is it too costly ?
    atom->cosPhase   = cos(2*M_PI*phase);
    atom->sinPhase   = sin(2*M_PI*phase);

    // Convert to HighResStft
    // The squared coefficient
    // Should take border type and stftHighRes->firstp/stftHighRes->lastp into effect ?
    SetAtomCoeff2HighRes(atom,stftHighRes,subStftReal,subStftPhase);

    // Store the computed data
    coeff2sHighRes[freqId/stftHighRes->fRate] = atom->coeff2;
  }
} 


/**************************************************************/
/*
 * Computes the values of the HighResStft stft associated with 
 * a given RealStft stft at the same windowSize and a RealStft stft
 * at a smaller windowSize (for High Time-Resolution) 
 * or (TODO) larger windowSize (for High Freq-Resolution). 
 * The signal associated with those RealStft stfts has changed
 * between the time-indexes [timeIdMin,timeIdMax].
 */
/**************************************************************/
void UpdateStftHighRes(STFT stftHighRes,
		       STFT stftReal,STFT stftPhase,
		       STFT subStftReal,STFT subStftPhase,
		       int timeIdMin,int timeIdMax)
{    
  int shiftMinAtom,shiftMaxAtom;
  int shiftMinSubAtom,shiftMaxSubAtom;
  int timeIdMinStft,timeIdMaxStft;
  //  int timeIdMinStft2,timeIdMaxStft2;
  //  char flagTwoLoops=NO;

  int timeId;
  
  /* Checking arguments */
  CheckStftHighResEnergy(stftHighRes,stftReal,stftPhase,subStftReal,subStftPhase);
  if((!INRANGE(0,timeIdMin,timeIdMax)) || 
     (!INRANGE(timeIdMin,timeIdMax,stftHighRes->signalSize-1)))
    Errorf("UpdateStftHighRes : bad time range [%d %d] not in [0 %d]",
	   timeIdMin,timeIdMax,stftHighRes->signalSize-1);
  
  
  // The support of a main atom at 'timeId' is
  // [timeId+shiftMinAtom,timeId+shiftMaxAtom]
  ComputeWindowSupport(stftHighRes->windowSize,stftHighRes->windowShape,0,
		     &shiftMinAtom,&shiftMaxAtom);
  // The support of a sub-atom at 'subTimeId' is
  // [subTimeId+shiftMinSubAtom,subTimeId+shiftMaxSubAtom]
  ComputeWindowSupport(subStftReal->windowSize,subStftReal->windowShape,0,
		     &shiftMinSubAtom,&shiftMaxSubAtom);

  // There are no boundary effects in the following region
  stftHighRes->firstp = shiftMaxSubAtom-shiftMinSubAtom-shiftMinAtom;
  stftHighRes->lastp  = stftHighRes->signalSize-1+shiftMinSubAtom-shiftMaxSubAtom-shiftMaxAtom;
  
  // The highres coefficient of a main atom can change
  // iff 
  //   timeId+shiftMaxAtom>=subTimeId+shiftMinSubAtom
  // where
  //   subTimeId+shiftMaxSubAtom>=timeIdMin
  // or
  //   timeId+shiftMinAtom<=subTimeId+shiftMaxSubAtom
  // where
  //   subTimeId+shiftMinSubAtom<=timeIdMax
  timeIdMinStft = timeIdMin-shiftMaxSubAtom
    +shiftMinSubAtom-shiftMaxAtom;
  timeIdMaxStft = timeIdMax-shiftMinSubAtom
    +shiftMaxSubAtom-shiftMinAtom;

  switch(stftHighRes->borderType) {
  case BorderPad0 :
    /* We consider only those on the  time-freq grid */
    QuantizeRangeLarge(timeIdMinStft,timeIdMaxStft,stftHighRes->tRate,
		       &timeIdMinStft,&timeIdMaxStft);
    timeIdMinStft = MAX(0,timeIdMinStft);
    timeIdMaxStft = MIN((stftHighRes->nFrames*stftHighRes->tRate)-stftHighRes->tRate,timeIdMaxStft);  
    break;
  default :
    Errorf("UpdateStftHighRes : border type %s not treated yet",BorderType2Name(stftHighRes->borderType));
  }

  // First time Loop 
  for (timeId = timeIdMinStft; 
       timeId <= timeIdMaxStft; 
       timeId += stftHighRes->tRate) {
    StftHighResSetFreq(stftHighRes,stftReal,stftPhase,subStftReal,subStftPhase,timeId);
  }
  stftHighRes->flagUpToDate = YES;
}

/* EOF */