StackOverflow2013

Note that there are some explanatory texts on larger screens.

plurals

POAlgorithm smbPitchShift (Pascal)
text
Body
copied!<p>I looked for a long time this algorithm in Pascal and not found, I found it only in C++, it was frustrating. Then I decided to translate the C++ code for Pascal, however there were some problems that I am not able to solve. it appeared an error message "Floating point overflow". I would like help to make this code work!</p> <pre><code>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)); </code></pre>

Querying!

Guidance

An individual column

Larger individual text columns get their own page to allow for proper reading.

SQuiL has stopped working due to an internal error.

If you are curious you may find further information in the browser console, which is accessible through the devtools (F12).

Reload