Clear Filters
Clear Filters

Multiple window calculation for FFT

3 views (last 30 days)
I`m analysing blood pressure using frequency domain. I have long recordings more than 20 min. I want to analyze it in 5 min-windows and get all the results of each window. However, when I get the results, I only get one value.
How can I get the all calculated values of band power and periodogram for each window??
% Analysis in 5 min, con 50% overlap
winsize = 300*fsr;
overlap=winsize/2;
j=1;
i = 1:overlap:(length(SAPr) - winsize);
% Define temporal series of each sample
SAPw = SAPr(i:(i + winsize)-1);
MAPw = MAPr(i:(i + winsize)-1);
DAPw = DAPr(i:(i + winsize)-1);
HRw = HRr(i:(i + winsize)-1);
% Dickey-Fuley test of stationarity
% if adftest(MAPw)==0 % 1 means stationarity
% continue;
% end
% Median values
meanSAP(j)= mean(SAPw);
meanMAP(j) = mean(MAPw);
meanDAP(j) = mean(DAPw);
meanHR(j) = mean(HRw);
% Detrend
SAPw = detrend(SAPw,10);
MAPw = detrend(MAPw,10);
DAPw= detrend(DAPw,10);
HRw = detrend(HRw,10);
%Standard deviation
stdSAP(j) = std(SAPw);
stdMAP(j) = std(MAPw);
stdDAP(j) = std(DAPw);
stdHR(j) = std(HRw);
% Periodogram
[PSAP,FS] = periodogram(SAPw,rectwin(length(SAPw)),length(SAPw),fsr);
[PMAP,FM] = periodogram(MAPw,rectwin(length(MAPw)),length(MAPw),fsr);
[PDAP,FD] = periodogram(DAPw,rectwin(length(DAPw)),length(DAPw),fsr);
[PHR,FH] = periodogram(HRw,rectwin(length(HRw)),length(HRw),fsr);
figure, plot(FS,PSAP); hold on; plot(FM,PMAP,'r');plot(FD,PDAP,'g');plot(FH,PHR,'k');
legend('SAP','DAP','MAP','HR');
% Band power
LFS(j) = bandpower(PSAP,FS,[0.033,0.15],'psd');
HFS(j) = bandpower(PSAP,FS,[0.151,0.4],"psd");
VLFS(j) = bandpower(PSAP,FS,[0,0.033],'psd');
TPS(j) = bandpower(PSAP,FS,[0,0.4],'psd');
LFrelS(j) = 100*(LFS/TPS);
HFrelS(j) = 100*(HFS/TPS);
LFHFS(j)= LFrelS(j)/HFrelS(j);
LFM(j) = bandpower(PMAP,FM,[0.033,0.15],'psd');
HFM(j) = bandpower(PMAP,FM,[0.151,0.4],"psd");
VLFM(j) = bandpower(PMAP,FM,[0,0.033],'psd');
TPM(j) = bandpower(PMAP,FM,[0,0.4],'psd');
LFrelM(j) = 100*(LFM/TPM);
HFrelM(j) = 100*(HFM/TPM);
LFHFM(j)= LFrelM(j)/HFrelM(j);
LFD(j) = bandpower(PDAP,FD,[0.033,0.15],'psd');
HFD(j) = bandpower(PDAP,FD,[0.151,0.4],"psd");
VLFD(j) = bandpower(PDAP,FD,[0,0.033],'psd');
TPD(j) = bandpower(PDAP,FD,[0,0.4],'psd');
LFrelD(j) = 100*(LFD/TPD);
HFrelD(j) = 100*(HFD/TPD);
LFHFD(j)= LFrelD(j)/HFrelD(j);
LFH(j) = bandpower(PHR,FH,[0.033,0.15],'psd');
HFH(j) = bandpower(PHR,FH,[0.151,0.4],"psd");
VLFH(j) = bandpower(PHR,FH,[0,0.033],'psd');
TPH(j) = bandpower(PHR,FH,[0,0.4],'psd');
LFrelH(j) = 100*(LFH/TPH);
HFrelH(j) = 100*(HFH/TPH);
LFHFH(j)= LFrelH(j)/HFrelH(j);
j=j+1
meanSAPu = mean(nonzeros(meanSAP));
meanMAPu = mean(nonzeros(meanMAP));
meanDAPu = mean(nonzeros(meanDAP));
meanHRu = mean(nonzeros(meanHR));
stdSAPu = mean(nonzeros(stdSAP));
stdMAPu = mean(nonzeros(stdMAP));
stdDAPu = mean(nonzeros(stdDAP));
stdHRu = mean(nonzeros(stdHR));
LFrelSu = mean(nonzeros(LFrelS));
LFrelMu = mean(nonzeros(LFrelM));
LFrelDu = mean(nonzeros(LFrelD));
LFrelHu = mean(nonzeros(LFrelH));
HFrelSu = mean(nonzeros(HFrelS));
HFrelMu = mean(nonzeros(HFrelM));
HFrelDu = mean(nonzeros(HFrelD));
HFrelHu = mean(nonzeros(HFrelH));
LFHFHu = mean(nonzeros(LFHFH));
LFHFSu = mean(nonzeros(LFHFS));
LFHFDu = mean(nonzeros(LFHFD));
LFHFMu = mean(nonzeros(LFHFM));
% Save results in an Excel file
results=table(meanSAP,stdSAP,LFrelS,HFrelS,meanMAP,stdMAP,LFrelM,HFrelM,meanDAP,stdDAP,LFrelD,HFrelD,meanHR,stdHR,LFrelH,HFrelH,LFHFS,LFHFM,LFHFD,LFHFH)
results2=splitvars(results)
resultadosfinal=rows2vars(results2)
writetable(resultadosfinal,"Prueba Matlab final.csv")

Answers (1)

Divit
Divit on 18 Dec 2023
Hi Tomas,
I understand that you are facing issues with window calculation and not getting the expected results.
It seems that the issue lies in how the loop is structured and how the indices are updated. The code provided is intended to process the signal in overlapping windows, but the loop only runs once because the indices are not updated correctly within the loop.
Here's the corrected code with a "for" loop that iterates over each segment of the signal:
% Analysis in 5 min, with 50% overlap
winsize = 300 * fsr;
overlap = winsize / 2;
% Preallocate arrays for efficiency
num_windows = floor((length(SAPr) - winsize) / overlap);
meanSAP = zeros(num_windows, 1);
meanMAP = zeros(num_windows, 1);
% ... (Do this for all other variables you're computing)
j = 1;
for i = 1:overlap:(length(SAPr) - winsize)
% Define temporal series of each sample
SAPw = SAPr(i:(i + winsize) - 1);
MAPw = MAPr(i:(i + winsize) - 1);
% ... (Do this for all other signals)
% Perform your analysis here
% ...
% Periodogram
[PSAP, FS] = periodogram(SAPw, rectwin(length(SAPw)), length(SAPw), fsr);
% ... (Do this for all other signals)
% Band power
LFS(j) = bandpower(PSAP, FS, [0.033, 0.15], 'psd');
% ... (Do this for all other bands and signals)
% Update the index for the results
j = j + 1;
end
% After the loop, you can calculate the mean of the non-zero elements if needed
meanSAPu = mean(nonzeros(meanSAP));
% ... (Do this for all other variables)
% Save results in an Excel file
results = table(meanSAP, stdSAP, LFrelS, HFrelS, meanMAP, stdMAP, LFrelM, HFrelM, meanDAP, stdDAP, LFrelD, HFrelD, meanHR, stdHR, LFrelH, HFrelH, LFHFS, LFHFM, LFHFD, LFHFH);
results2 = splitvars(results);
resultadosfinal = rows2vars(results2);
writetable(resultadosfinal, "Prueba Matlab final.csv");
Here's what was changed:
  1. The loop is now a "for" loop that iterates over the starting indices of each window, with the correct step size to account for the overlap.
  2. The "j" index is initialized outside the loop and incremented inside the loop to store the results for each window.
  3. The preallocation of result arrays to improve performance and ensure that they are of the correct size.
  4. The calculation of mean values for non-zero elements and the saving of results to an Excel file are done after the loop, ensuring that all windowed results are included.
With these changes, the code should now process each 5-minute window of your signal and store the results in the corresponding arrays. After the loop, you can calculate the mean across all windows if needed and then save the results to a file.
I hope it helps!

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!