2011-08-31 15 views
5

Busqué durante mucho tiempo este algoritmo en Pascal y no lo encontré, lo encontré solo en C++, fue frustrante. Luego decidí traducir el código C++ para Pascal, sin embargo, hubo algunos problemas que no puedo resolver. apareció un mensaje de error "Desbordamiento de punto flotante". ¡Me gustaría ayuda para hacer que este código funcione!Algoritmo smbPitchShift (Pascal)

var 
    WFX: pWaveFormatEx; 

    {** Algoritimo Pitch Shift **} 
    gInFIFO, gOutFIFO, gLastPhase, gSumPhase, gOutputAccum: Array Of Extended; 
    gAnaMagn, gAnaFreq, gSynFreq, gSynMagn, gFFTworksp: Array Of Extended; 

Const 
    MAX_FRAME_LENGTH = 8192; 

implementation 

{$R *.dfm} 

procedure smbFft(fftBuffer: PExtended; fftFrameSize, sign: Integer); 
var 
    p1, p2, p1r, p1i, p2r, p2i: PExtended; 
    wr, wi, arg, temp: EXTENDED; 
    tr, ti, ur, ui: EXTENDED; 
    i, bitm, j, le, le2, k: Integer; 
begin 
    i:= 2; 
    WHILE (i < 2*fftFrameSize-2) DO            //for (i = 2; i < 2*fftFrameSize-2; i += 2) { 
    BEGIN 
     bitm:= 2; 
     j:= 0; 
     WHILE (bitm < (2 * fftFrameSize)) DO          //for (bitm = 2, j = 0; bitm < 2*fftFrameSize; bitm <<= 1) { 
     BEGIN 
     if ((i and bitm) <> 0) then           //if (i & bitm) j++; 
      inc(j); 
      // 
      j:= j shl 1;               //j <<= 1; 
      bitm:= bitm shl 1;              //bitm <<= 1 
    END; 
     // 
     if (i < j) then 
     begin 
     p1:= fftBuffer;              //^ 
      Inc(p1, i);               //p1 = fftBuffer+i; 
      p2:= fftBuffer;              //^ 
      Inc(p2, j);               //p2 = fftBuffer+j; 
     temp:= p1^;               //temp = *p1; 
      inc(p1, 1);               //*(p1++) 
      p1:= p2;                //p1 = *p2; 
      inc(p2, 1);               //*(p2++) 
      p2^:= temp;               //p2 = temp; 
      temp:= p1^;               //temp = *p1; 
     p1:= p2;                //*p1 = *p2; 
      p2^:= temp;               //*p2 = temp; 
    end; 
     INC(I, 2); 
    END; 
    // 
    le:= 2; 
    k:= 0; 
    WHILE (k < (ln(fftFrameSize)/ln(2.0)+0.5)) DO         //for (k = 0, le = 2; k < (long)(log(fftFrameSize)/log(2.)+.5); k++) { 
    BEGIN 
     le:= le shl 1;                //le <<= 1; 
     le2:= le shr 1;               //le2 = le>>1; 
     ur:= 1.0;                 //ur = 1.0; 
     ui:= 0.0;                 //ui = 0.0; 
     arg:= PI/(le2 shr 1);             //arg = M_PI/(le2>>1); 
     wr:= cos(arg);                //wr = cos(arg); 
     wi:= sign * sin(arg);              //wi = sign*sin(arg); 
     j:=0; 
     WHILE (j < le2) DO               //for (j = 0; j < le2; j += 2) { 
     BEGIN 
      p1r:= fftBuffer;              //^ 
      INC(p1r, j);               //p1r = fftBuffer+j; 
      p1i:= p1r;                //^ 
      INC(p1i, 1);               //p1i = p1r+1; 
      p2r:= p1r;                //^ 
      INC(p2r, le2);               //p2r = p1r+le2; 
      p2i:= p2r;                //^ 
      INC(p2i, 1);               //p2i = p2r+1; 
      i:= j; 
     WHILE (i < 2*fftFrameSize) DO           //for (i = j; i < 2*fftFrameSize; i += le) { 
      BEGIN 
      tr:= p2r^ * ur - p2i^ * ui;          //tr = *p2r * ur - *p2i * ui; 
      ti:= p2r^ * ui + p2i^ * ur;          //ti = *p2r * ui + *p2i * ur; 
      p2r^:= p1r^ - tr;             //*p2r = *p1r - tr; 
       p2i^:= p1i^ - ti;             //*p2i = *p1i - ti; 
      p1r^:= p1r^ + tr;             //*p1r += tr; 
       p1i^:= p1i^ + ti;             //*p1i += ti; 
      INC(p1r, le);              //p1r += le; 
       INC(p1i, le);              //p1i += le; 
       INC(p2r, le);              //p2r += le; 
       INC(p2i, le);              //p2i += le; 
       INC(i, le); 
     END; 
      // 
      tr:= ur * wr - ui * wi;            //tr = ur*wr - ui*wi; 
     ui:= ur * wi + ui * wr;            //ui = ur*wi + ui*wr; 
      ur:= tr;                //ur = tr; 
      INC(J, 2); 
      END; 
     inc(k); 
    END; 
end; 

Procedure smbPitchShift(pitchShift: Double; numSampsToProcess, fftFrameSize, osamp, sampleRate: Integer; indata, outdata: PExtended); 

    function atan2 (y, x : Extended) : Extended; Assembler; 
    asm 
    fld [y] 
    fld [x] 
    fpatan 
    end; 
var magn, phase, tmp, window, xreal, imag: Extended; 
    freqPerBin, expct, CC: Extended; 
    i, k, qpd, index, inFifoLatency, stepSize, fftFrameSize2: Integer; 
    gRover: Integer; 
    TmpData: PExtended; 
begin 
    gRover:= 0; 
    {* set up some handy variables *} 
    fftFrameSize2:= Round(fftFrameSize/2);          //fftFrameSize2 = fftFrameSize/2; 
    stepSize:= Round(fftFrameSize/osamp);          //stepSize = fftFrameSize/osamp; 
    freqPerBin:= sampleRate/fftFrameSize;          //freqPerBin = sampleRate/(double)fftFrameSize; 
    expct:= 2.0 * PI * stepSize/fftFrameSize;         //expct = 2.*M_PI*(double)stepSize/(double)fftFrameSize; 
    inFifoLatency:= fftFrameSize - stepSize;          //inFifoLatency = fftFrameSize-stepSize; 
    if (gRover = 0) then gRover:= inFifoLatency;         //if (gRover == false) gRover = inFifoLatency; 
    // 
    {* main processing loop *} 
    for i:=0 to numSampsToProcess-1 do            //for (i = 0; i < numSampsToProcess; i++){ 
    begin 
     {* As long as we have not yet collected enough data just read in *} 
     TmpData:= indata;               //^ 
     inc(TmpData, i);               // [i] 
     gInFIFO[gRover]:= TmpData^;            //gInFIFO[gRover] = indata[i]; 
     TmpData:= outdata;               //^ 
     inc(TmpData, i);               // [i] 
     TmpData^:= gOutFIFO[gRover - inFifoLatency];        //outdata[i] = gOutFIFO[gRover-inFifoLatency]; 
     Inc(gRover);                //gRover++; 
     {* now we have enough data for processing *} 
     if (gRover >= fftFrameSize) then           //if (gRover >= fftFrameSize) { 
      begin 
      gRover:= inFifoLatency;            //gRover = inFifoLatency; 
      {* do windowing and re,im interleave *} 
      for k:=0 to fftFrameSize-1 do          //for (k = 0; k < fftFrameSize;k+ 
       begin 
       window:= -0.5 * Cos(2.0 * PI * k/fftFrameSize) + 0.5;   //window = -.5*cos(2.*M_PI*(double)k/(double)fftFrameSize)+.5; 
       gFFTworksp[2 * k]:= gInFIFO[k] * window;       //gFFTworksp[2*k] = gInFIFO[k] * window; 
       gFFTworksp[2 * k + 1]:= 0.0;          //gFFTworksp[2 * k + 1]:= 0.0F; 
       end; 
      {****************** ANALYSIS ********************} 
      {* do transform *} 
      SmbFft(Ptr(DWORD(gFFTworksp)), fftFrameSize, -1);     //smbFft(gFFTworksp, fftFrameSize, -1); 
      {* this is the analysis step *} 
      for k:= 0 to fftFrameSize2 do          //for (k = 0; k <= fftFrameSize2; k++) { 
       begin 
       {* de-interlace FFT buffer *} 
       xreal:= gFFTworksp[2 * k];          //real = gFFTworksp[2*k]; 
       imag:= gFFTworksp[2 * k + 1];         //imag = gFFTworksp[2*k+1]; 
       {* compute magnitude and phase *} 
       magn:= 2.0 * Sqrt(xreal * xreal + imag * imag);     //magn = 2.*sqrt(real*real + imag*imag); 
       phase:= Atan2(imag, xreal);          //phase = atan2(imag,real); 
       {* compute phase difference *} 
       tmp:= phase - gLastPhase[k];          //tmp = phase - gLastPhase[k]; 
       gLastPhase[k]:= phase;           //gLastPhase[k] = phase; 
       {* subtract expected phase difference *} 
       tmp:= tmp - k * expct;           //tmp -= (double)k*expct; 
       {* map delta phase into +/- Pi interval *} 
       qpd:= Round(tmp/PI);           //qpd = tmp/M_PI; 
       if (qpd >= 0) then 
        qpd:= qpd + qpd and 1           // if (qpd >= 0) qpd += qpd&1; 
       else 
        qpd:= qpd - qpd and 1;           // else qpd -= qpd&1; 
       // 
       tmp:= tmp - (PI * qpd);           //tmp -= M_PI*(double)qpd; 
       {* get deviation from bin frequency from the +/- Pi interval *} 
       tmp:= osamp * tmp/(2.0 * PI);         //tmp = osamp*tmp/(2.*M_PI); 
       {* compute the k-th partials' true frequency *} 
       tmp:= k * freqPerBin + tmp * freqPerBin;       //tmp = (double)k*freqPerBin + tmp*freqPerBin; 
       {* store magnitude and true frequency in analysis arrays *} 
       gAnaMagn[k]:= magn;            //gAnaMagn[k] = magn; 
       gAnaFreq[k]:= tmp;            //gAnaFreq[k] = tmp; 
       end; 
      {****************** PROCESSING ********************} 
      {* this does the actual pitch shifting *} 
      for k:=0 to fftFrameSize2 do           //for (k = 0; k <= fftFrameSize2; k++) { 
       begin 
       index:= Round(k * pitchShift);         //index = (long)(k*pitchShift); 
       if (index <= fftFrameSize2) then         //if (index <= fftFrameSize2) { 
        begin 
        IF K >= LENGTH(gSynFreq) THEN 
         SetLength(gSynFreq , LENGTH(gSynFreq)+1);     //memset(gSynFreq, 0, fftFrameSize*sizeof(float)); 
        IF K >= LENGTH(gSynMagn) THEN 
         SetLength(gSynMagn , LENGTH(gSynMagn)+1);     //memset(gSynMagn, 0, fftFrameSize*sizeof(float)); 
        // 
        gSynMagn[index]:= gSynMagn[index] + gAnaMagn[k];    //gSynMagn[index] += gAnaMagn[k]; 
        gSynFreq[index]:= gAnaFreq[k] * pitchShift;     //gSynFreq[index] = gAnaFreq[k] * pitchShift; 
        end; 
       end; 
      {****************** SYNTHESIS ********************} 
      {* this is the synthesis step *} 
      for k:=0 to fftFrameSize2 do           //for (k = 0; k <= fftFrameSize2; k++) { 
       begin 
       {* get magnitude and true frequency from synthesis arrays *} 
       magn:= gSynMagn[k];            // magn = gSynMagn[k]; 
       tmp:= gSynFreq[k];            //tmp = gSynFreq[k] 
       {* subtract bin mid frequency *} 
       tmp:= tmp - (k * freqPerBin);         //tmp -= (double)k*freqPerBin; 
       {* get bin deviation from freq deviation *} 
       tmp:= tmp/freqPerBin;           //tmp /= freqPerBin; 
       {* take osamp into account *} 
       tmp:= 2.0 * PI * tmp/osamp;         //tmp = 2.*M_PI*tmp/osamp; 
       {* add the overlap phase advance back in *} 
       tmp:= tmp + (k * expct);           //tmp += (double)k*expct; 
       {* accumulate delta phase to get bin phase *} 
       gSumPhase[k]:= gSumPhase[k] + tmp;        //gSumPhase[k] += tmp; 
       phase:= gSumPhase[k];           //phase = gSumPhase[k]; 
       {* get real and imag part and re-interleave *} 
       gFFTworksp[2 * k]:= (magn * Cos(phase));       //gFFTworksp[2*k] = magn*cos(phase); 
       gFFTworksp[2 * k + 1]:= (magn * Sin(phase));      //gFFTworksp[2*k+1] = magn*sin(phase); 
       end; 
      {* zero negative frequencies *} 
      k:= fftFrameSize + 2; 
      WHILE (k < 2 * fftFrameSize) DO          //for (k = fftFrameSize+2; k < 2*fftFrameSize; k++) 
       BEGIN 
       gFFTworksp[k]:= 0.0;            //gFFTworksp[k] = 0.0F; 
       inc(k); 
       END; 
      {* do inverse transform *} 
      SmbFft(Ptr(DWORD(gFFTworksp)), fftFrameSize, 1);      //smbFft(gFFTworksp, fftFrameSize, 1); 
      {* do windowing and add to output accumulator *} 
      for k:=0 to fftFrameSize-1 do          // for(k=0; k < fftFrameSize; k++) { 
       begin 
       window:= -0.5 * Cos(2.0 * PI * k/fftFrameSize) + 0.5;   //window = -.5*cos(2.*M_PI*(double)k/(double)fftFrameSize)+.5; 
       gOutputAccum[k]:= gOutputAccum[k] + (2.0 * window * gFFTworksp[2 * k]/(fftFrameSize2 * osamp)); 
       end;                //gOutputAccum[k] += 2.*window*gFFTworksp[2*k]/(fftFrameSize2*osamp); 
      // 
      for k:=0 to stepSize-1 do gOutFIFO[k]:= gOutputAccum[k];    //for (k = 0; k < stepSize; k++) gOutFIFO[k] = gOutputAccum[k]; 
      {* shift accumulator *} 
      // 
      TmpData:= PTR(DWORD(gOutputAccum));         //^ 
      Inc(TmpData, StepSize);            //gOutputAccum + stepSize 
      MoveMemory(TmpData, PTR(DWORD(gOutputAccum)), fftFrameSize * sizeof(Extended)); 
                      //memmove(gOutputAccum, gOutputAccum + stepSize, fftFrameSize * sizeof(float)); 
      // 
      {* move input FIFO *} 
      for k:=0 to inFifoLatency-1 do          //for (k = 0; k < inFifoLatency; k++) 
       gInFIFO[k]:= gInFIFO[k + stepSize];        //gInFIFO[k] = gInFIFO[k+stepSize]; 
      end; 
    end; 
end; 

procedure TWavAnalize.FormCreate(Sender: TObject); 
begin 
    {** algoritimo pitchshift **} 
    SetLength(gInFIFO ,MAX_FRAME_LENGTH); 
    SetLength(gOutFIFO ,MAX_FRAME_LENGTH); 
    SetLength(gSynFreq ,MAX_FRAME_LENGTH); 
    SetLength(gSynMagn ,MAX_FRAME_LENGTH); 
    SetLength(gAnaFreq ,MAX_FRAME_LENGTH); 
    SetLength(gAnaMagn ,MAX_FRAME_LENGTH); 
    SetLength(gFFTworksp ,2 * MAX_FRAME_LENGTH); 
    SetLength(gLastPhase , Round(MAX_FRAME_LENGTH/2) + 1); 
    SetLength(gSumPhase , Round(MAX_FRAME_LENGTH/2) + 1); 
    SetLength(gOutputAccum , 2 * MAX_FRAME_LENGTH); 
    {** algoritimo pitchshift **} 
end; 

procedure TWavAnalize.Button8Click(Sender: TObject); 
VAR T: TMEMORYSTREAM; 
    DSize, DataOffset, i: cARDINAL; 
    AIN, AOUT: ARRAY OF EXTENDED; 
begin 
    T:= TMEMORYSTREAM.CREATE; 
    T.LoadFromFile(PATH); 
    GetStreamWaveAudioInfo(T, WFX, DSize, DataOffset); 
    T.Position:= DataOffset; 
    SETLENGTH(AIN, DSIZE); 
    SETLENGTH(AOUT, DSIZE); 
    T.READ(AIN[0], DSIZE); 
    smbPitchShift(0.5, DSize, 2048, 10, WFX.nSamplesPerSec, Ptr(DWORD(AIN)), Ptr(DWORD(AOUT))); 
    T.Clear; 
    T.WRITE(AOUT[0], LENGTH(AOUT)); 
+1

Bienvenido a StackOverflow. Por favor, edite su pregunta para formatear el código correctamente. Utilice el botón con la imagen '{}' en la barra de herramientas o Ctrl + K para formatearla una vez que esté pegada. Puede obtener una vista previa de su publicación debajo de donde se ingresa simplemente mirando esa área de la ventana del navegador. Iba a tratar de arreglarlo para ti, pero es demasiado complicado (y falta parte al final) para que lo haga. Si las personas no pueden leerlo, no pueden ayudarte. –

+3

FWIW, en C y C++ en Win32 y Win64, *** float *** significa *** Single ***, no Extended. No es necesario usar Extended en el código. –

+1

¿Por qué no usa la DLL de SoundTouch como le dije en una de sus otras preguntas sobre el cambio de tono? –

Respuesta

3

Tenga en cuenta los HINTS que Delphi genera al compilar el código.

Por ejemplo me sale: "[DCC Sugerencia] Unit1.pas (65): Valor asignado a H2077 'P1' nunca utilizado" en esta línea:

p1:= p2;                //*p1 = *p2; 

Debido a que en realidad debería ser:

p1^ := p2^; 

Además, si se agrega esta línea en la parte superior de la unidad:

{$POINTERMATH ON} 

que puede hacer Pointe C-estilo r aritmética por lo que no tiene que usar la solución TmpData. Así que esto:

TmpData:= indata;             
inc(TmpData, i);              
gInFIFO[gRover]:= TmpData^;  

se puede simplificar en esto:

gInFIFO[gRover]:= InData[i]; 
6

OK, por lo general no hacen esto, pero también tengo un interés en tener una versión de Delphi del código, por lo que lo tradujo Pruebe my translation y vea si eso funciona para usted.

FWIW, también eche un vistazo a Dirac3LE library del mismo autor. Esa es una biblioteca mucho más profesional (PSOLA, no WSOLA), disponible para Windows, Linux, Mac, iPhone, etc. Acabo de probar la versión para Mac y suena bien.

+0

Buen trabajo, así debería ser una traducción. –

+0

@Ville: Gracias. No es mi primer.

+0

+1 Bien, jugaré con esa noche –