function [ meanMSE ] = genMeanMSEResultAWGN( filtParams, FFTblocks, filtFunc, numNoiseRealizations,...
                                                     X,Y,noiseMag, imgParams, newLfInv  )
%UNTITLED2 Summary of this function goes here
%   Detailed explanation goes here


% tic();

[numRows,numCols] = size(newLfInv);

[numX,numY]=size(X);

fftNoiseMult = sqrt(numX*numY);

numChannels = size(FFTblocks,3);



% for k = 1:numChannels
%     F{k} = planckTaper2D(filtParams(1),filtParams(2),X,Y);    
F = planckTaper2D(filtParams(1),filtParams(2),X,Y);    
% end
% tic();
% for l=1:numNoiseRealizations    
%     for k = 1:numRows
%         tempFFTImg = zeros(numX,numY);
%         for j = 1:numCols    
%             if noiseMag > 0
%                 randNoise = imgParams.numPElectrons.*noiseMag.*fftNoiseMult.*randn(numX,numY);          
%             else
%                 randNoise = 0;
%             end 
%             tempFFTImg = tempFFTImg + newLfInv(k,j).*F.*(FFTblocks(:,:,j) + randNoise);                
%         end                
%         fftImg = tempFFTImg./(imgParams.imgTransFactor(k)*imgParams.imgMult)-imgParams.Zfft(:,:,k);
%         MSEtemp(l,k) = sum(sum(abs(fftImg).^2))./(numX*numY);
%     end    
% end
for l=1:numNoiseRealizations    
    if noiseMag > 0
        randNoise = imgParams.numPElectrons.*noiseMag.*fftNoiseMult.*randn(numX,numY,numRows);        
    else
        randNoise = 0;
    end
    FFTblocksNoise = FFTblocks + randNoise;
    for k = 1:numCols        
        %filtFFTBlocks = bsxfun(@times,FFTblocksNoise,F);
       filtFFTBlocks(:,:,k) = F.*FFTblocksNoise(:,:,k);        
    end
%     filtFFTBlocks = permute(filtFFTBlocks,[3,1,2,4]);
    %mmx cannot handle imaginary units
%     tempFFTImgTot = mmxc('mult',newLfInv,filtFFTBlocks);
%     tempFFTImgTot = zeros(numRows,numX,numY,numNoiseRealizations);
    tempFFTImgTot = zeros(numX,numY,numRows);
    for k = 1:numRows        
%         tic();
        for j = 1:numCols                
            tempFFTImgTot(:,:,k) = tempFFTImgTot(:,:,k) + newLfInv(k,j).*filtFFTBlocks(:,:,j);                
        end   
    end    
%     tempFFTImgTot = bsxfun(@times,tempFFTImgTot,(1./imgParams.imgTransFactor).');
    tempFFTImgTot = tempFFTImgTot./(imgParams.imgMult);
    fftImgTot = tempFFTImgTot-imgParams.Zfft;
    MSEtemp(l,:) = sum(sum(abs(fftImgTot).^2,1),2)./(numX*numY);
end
%toc()
% geometric mean between filter types to normalize differences in the data 
% between filter types for the cost function.
meanMSE = geomean(mean(MSEtemp,1)); 
% toc()

end

