function p= HHTobradaMNF (ctrlPrikaza,signal,IMFcijeli,kut,Fs,kutLMA,L,R,LMA)
% function p= HHTobradaMNF (ctrlPrikaza,signal,IMFcijeli,kut,Fs,kutLMA,L,R,LMA)
%
% Obrada za L kao parametar, uz R=L
% IMFs JESU prethodno odreðeni

MNF= []; 
tMNF= [];


for i=1:R:size(IMFcijeli,2)-L
    
    imf= [];
    
    imf= IMFcijeli(:,i:i+L-1); % bolji rezultati za novu EMD funkciju(FLANDRIN)
    tMNF{end+1}= (i+i+L-1)/2;
    
 
    
    %-->###################################################################
    % raèunanje HILBERTOVOG spektra
    % A=amplituda, f=frekvencija(normirana)
    %
    [A,f,tt]= hhspectrum(imf);
    %[A,f,tt]= hhspectrum2(imf);
    f= f*Fs; % množimo s Fs jer je frekvencija normirana
    %###################################################################<--
    
    
    
    %-->###################################################################
    % izraèunavanje težinske srednje frekvencije za svaki IMF
    % rezultat pohranjen u MIF
    MIF= [];

    for k = 1:size(A,1)

        MIF{end+1}= sum(f(k,:).*A(k,:).*A(k,:))/sum(A(k,:).*A(k,:));   
        
    end;
    %###################################################################<--
    
    
    
    %-->###################################################################
    % raèunanje frekvencije medijana vremenskog odsjeèka na naèin da
    % frekvenciju medijana svakog IMF težinski usrednjimo dvostrukom normom 
    % svakog pojedinog IMFa.
    % rezultat pohranjen u MNF.
    An= [];

    for k = 1:size(A,1)
    
        %A= Amp{k};
        An{end+1}= sum(A(k,:).*A(k,:))^0.5;
               
    end;
    
    MIF= cell2mat(MIF);
    An= cell2mat(An);
    
    % izraèunata frekvencija medijana za pojedini vremenski odsjeèak
    MNF{end+1}= sum(MIF.*An)/sum(An); 
    %###################################################################<--
    
    
%     clear imf;
%     clear MIF;
%     clear An;
%         
%     clear A;
%     clear imf1;
    
       
end   



%% Izraèunavanje REGRESIJSKOG pravca

% izglaðivnaje izraèunatih frekvencija sa dvostrukim filtriranjem.
HHTsmMNF= smoothMF(cell2mat(MNF),LMA);

% odreðivanje iznosa i vremenskih trenutaka lokalnih maksimuma
[lokMAX, tMAX]= lokMaximum(HHTsmMNF,cell2mat(tMNF),kutLMA,kut);

% linearna interpolacija u smislu najmanje kvadratne pogreške -LMSE
% pREG je regresijski pravac izraèunat za vremensku skalu u [ms].
% za skalu u [sekundama] treba množiti p(1) sa 1000.
% pREG= p(1)*x + p(2)
p= polyfit(tMAX,lokMAX,1);
pREG= p(1)*[1:length(signal)]+ p(2);


%% PRIKAZ dobivenih velièina

if ctrlPrikaza==1
    
    figure,
    subplot(4,1,1); plot([1:length(signal)]/1000,signal/1000);
        xlabel('vrijeme [s]'); ylabel('[mV]'); 
        title(['EMG signal: L= ',num2str(L),', R= ',num2str(R)]);
    subplot(4,1,2); plot(cell2mat(tMNF)/1000,cell2mat(MNF)); 
        xlabel('vrijeme [s]'); ylabel('[Hz]'); 
        title('Srednja frekvencija - MNF');
    subplot(4,1,3); plot(cell2mat(tMNF)/1000,HHTsmMNF); hold on;
        plot(tMAX/1000,lokMAX,'rx','MarkerSize',10); 
        xlabel('vrijeme [s]'); ylabel('[Hz]'); 
        title('Izglaðena srednja frekvencija - MNF'); 
        hold off;
    subplot(4,1,4); plot([1:length(signal)]/1000,pREG,'k'); hold on; 
        plot(tMAX/1000,lokMAX,'rx','MarkerSize',10); 
        xlabel('vrijeme [s]'); ylabel('[Hz]'); 
        title('Regresijski pravac i lokalni maksimumi');
        str1= strcat('y= ',num2str(60*1000*p(1)),'[Hz/min]*x + ',num2str(p(2)));
        str2= 'lok. maksimumi';
        legend(str1,str2);
        hold off;

end