I searched for this algorithm in Pascal for a long time and did not find it, 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 could not solve. I got the error "Floating point overflow". I would like to help make this code work!
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
BEGIN
bitm:= 2;
j:= 0;
WHILE (bitm < (2 * fftFrameSize)) DO
BEGIN
if ((i and bitm) <> 0) then
inc(j);
j:= j shl 1;
bitm:= bitm shl 1;
END;
if (i < j) then
begin
p1:= fftBuffer;
Inc(p1, i);
p2:= fftBuffer;
Inc(p2, j);
temp:= p1^;
inc(p1, 1);
p1:= p2;
inc(p2, 1);
p2^:= temp;
temp:= p1^;
p1:= p2;
p2^:= temp;
end;
INC(I, 2);
END;
le:= 2;
k:= 0;
WHILE (k < (ln(fftFrameSize)/ln(2.0)+0.5)) DO
BEGIN
le:= le shl 1;
le2:= le shr 1;
ur:= 1.0;
ui:= 0.0;
arg:= PI / (le2 shr 1);
wr:= cos(arg);
wi:= sign * sin(arg);
j:=0;
WHILE (j < le2) DO
BEGIN
p1r:= fftBuffer;
INC(p1r, j);
p1i:= p1r;
INC(p1i, 1);
p2r:= p1r;
INC(p2r, le2);
p2i:= p2r;
INC(p2i, 1);
i:= j;
WHILE (i < 2*fftFrameSize) DO
BEGIN
tr:= p2r^ * ur - p2i^ * ui;
ti:= p2r^ * ui + p2i^ * ur;
p2r^:= p1r^ - tr;
p2i^:= p1i^ - ti;
p1r^:= p1r^ + tr;
p1i^:= p1i^ + ti;
INC(p1r, le);
INC(p1i, le);
INC(p2r, le);
INC(p2i, le);
INC(i, le);
END;
tr:= ur * wr - ui * wi;
ui:= ur * wi + ui * wr;
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);
stepSize:= Round(fftFrameSize / osamp);
freqPerBin:= sampleRate / fftFrameSize;
expct:= 2.0 * PI * stepSize / fftFrameSize;
inFifoLatency:= fftFrameSize - stepSize;
if (gRover = 0) then gRover:= inFifoLatency;
{* main processing loop *}
for i:=0 to numSampsToProcess-1 do
begin
{* As long as we have not yet collected enough data just read in *}
TmpData:= indata;
inc(TmpData, i);
gInFIFO[gRover]:= TmpData^;
TmpData:= outdata;
inc(TmpData, i);
TmpData^:= gOutFIFO[gRover - inFifoLatency];
Inc(gRover);
{* now we have enough data for processing *}
if (gRover >= fftFrameSize) then
begin
gRover:= inFifoLatency;
{* do windowing and re,im interleave *}
for k:=0 to fftFrameSize-1 do
begin
window:= -0.5 * Cos(2.0 * PI * k / fftFrameSize) + 0.5;
gFFTworksp[2 * k]:= gInFIFO[k] * window;
gFFTworksp[2 * k + 1]:= 0.0;
end;
{****************** ANALYSIS ********************}
{* do transform *}
SmbFft(Ptr(DWORD(gFFTworksp)), fftFrameSize, -1);
{* this is the analysis step *}
for k:= 0 to fftFrameSize2 do
begin
{* de-interlace FFT buffer *}
xreal:= gFFTworksp[2 * k];
imag:= gFFTworksp[2 * k + 1];
{* compute magnitude and phase *}
magn:= 2.0 * Sqrt(xreal * xreal + imag * imag);
phase:= Atan2(imag, xreal);
{* compute phase difference *}
tmp:= phase - gLastPhase[k];
gLastPhase[k]:= phase;
{* subtract expected phase difference *}
tmp:= tmp - k * expct;
{* map delta phase into +/- Pi interval *}
qpd:= Round(tmp / PI);
if (qpd >= 0) then
qpd:= qpd + qpd and 1
else
qpd:= qpd - qpd and 1;
tmp:= tmp - (PI * qpd);
{* get deviation from bin frequency from the +/- Pi interval *}
tmp:= osamp * tmp / (2.0 * 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));