เร่งการคำนวณ FFT แบบกระจัดกระจาย

ฉันหวังว่าจะมีคนตรวจสอบโค้ดของฉันด้านล่างและเสนอคำแนะนำวิธีเร่งความเร็วส่วนระหว่าง tic และ toc ฟังก์ชันด้านล่างพยายามดำเนินการ IFFT เร็วกว่าฟังก์ชันในตัวของ Matlab เนื่องจาก (1) ถังขยะค่าสัมประสิทธิ์ fft เกือบทั้งหมดเป็นศูนย์ (เช่น ถังขยะ 10 ถึง 1000 จาก 10M ถึง 300M ถังขยะไม่เป็นศูนย์) และ (2 ) มีเพียงผลลัพธ์ IFFT ส่วนตรงกลางเท่านั้นที่ยังคงอยู่ (ผลลัพธ์แรกและสามสุดท้ายจะถูกละทิ้ง - ดังนั้นจึงไม่จำเป็นต้องคำนวณตั้งแต่แรก)

ตัวแปรอินพุตคือ:

fftcoef = complex fft-coef 1D array (10 to 1000 pts long)
bins = index of fft coefficients corresponding to fftcoef (10 to 1000 pts long)
DATAn = # of pts in data before zero padding and fft (in range of 10M to 260M)
FFTn = DATAn + # of pts used to zero pad before taking fft (in range of 16M to 268M) (e.g. FFTn = 2^nextpow2(DATAn))

ในปัจจุบัน โค้ดนี้ใช้เวลาสองสามลำดับของขนาดนานกว่าวิธีฟังก์ชัน ifft ของ Matlab ซึ่งจะคำนวณสเปกตรัมทั้งหมด จากนั้นจึงละทิ้ง 2/3 ของมัน ตัวอย่างเช่น หากข้อมูลอินพุตสำหรับ fftcoef และ bins เป็นอาร์เรย์ 9x1 (นั่นคือ เฉพาะค่าสัมประสิทธิ์ fft ที่ซับซ้อน 9 ต่อแถบด้านข้างเท่านั้น; 18 pts เมื่อพิจารณาแถบด้านข้างทั้งสอง) และ DATAn=32781534, FFTn=33554432 (เช่น 2^25) ดังนั้นวิธี ifft จะใช้เวลา 1.6 วินาที ในขณะที่ วนซ้ำด้านล่างใช้เวลามากกว่า 700 วินาที

ฉันหลีกเลี่ยงการใช้เมทริกซ์เพื่อ vectorize nn loop เนื่องจากบางครั้งขนาดอาเรย์สำหรับ fftcoef และ bins อาจยาวได้ถึง 1000 pts และเมทริกซ์ 260Mx1K จะมีขนาดใหญ่เกินไปสำหรับหน่วยความจำ เว้นแต่ว่ามันจะถูกทำลายด้วยวิธีใดวิธีหนึ่ง

คำแนะนำใด ๆ ที่ชื่นชมมาก! ขอบคุณล่วงหน้า.

function fn_fft_v1p0(fftcoef, bins, DATAn, FFTn)

fftcoef = [fftcoef; (conj(flipud(fftcoef)))];     % fft coefficients
bins = [bins; (FFTn - flipud(bins) +2)];          % corresponding fft indices for fftcoef array

ttrend = zeros( (round(2*DATAn/3) - round(DATAn/3) + 1), 1); % preallocate

start = round(DATAn/3)-1;

tic;
for nn = start+1 : round(2*DATAn/3)  % loop over desired time indices
  % sum over all fft indices having non-zero coefficients
  arg = 2*pi*(bins-1)*(nn-1)/FFTn;
  ttrend(nn-start) = sum( fftcoef.*( cos(arg) + 1j*sin(arg)); 
end
toc;

end