Функция нахождения среднеквадратического отклонения для оценки качества изображения

Image=imread('J:\diplom\512.bmp'); % читаем изображение
% size(I);
% whos I; 
%I=I(:,:,1);
Image1=rgb2gray(Image); %преобразуем цветное изображение в полутоновое и получаем его двумерный массив
M=im2double(Image1); %получаем трехмерный массив B содержащий изображение
%Z=im2unit8(M);
figure, imshow(M);title('Исходное изображение');
H=hadamard(512);
F=H*M; %прямое преобразование Адамара
W=H*F; %обратное преобразование ч.1
Y=W/512; %обратное преобразование с делением на N ч.2
X=im2uint8(Y);
figure, imshow(X);title('Изображение после обратного преобразования');

Function s = std2(a)
%STD2 Standard deviation of matrix elements.
%   B = STD2(A) computes the standard deviation of the values in
%   A.
%
%   Class Support
%   -------------
%   A can be numeric or logical. B is a scalar of class double.
%
%   Example
%   -------
%       I = imread('liftingbody.png');
%       val = std2(I)
%
%   See also CORR2, MEAN2, MEAN, STD.
 
%   Copyright 1992-2010 The MathWorks, Inc.
%   $Revision: 5.19.4.6 $  $Date: 2010/09/13 16:14:04 $
 
% validate that our input is valid for the IMHIST optimization
fast_data_type = isa(a,'logical') || isa(a,'int8') || isa(a,'uint8') || ...
    isa(a,'uint16') || isa(a,'int16');
 
% only use IMHIST for images of sufficient size
big_enough = numel(a) > 300000;
 
if fast_data_type && isequal(ndims(a),2) && ~issparse(a) && big_enough
    
    % compute histogram
    if islogical(a)
        num_bins = 2;
    else
        data_type = class(a);
        num_bins = double(intmax(data_type)) - double(intmin(data_type)) + 1;
    end
    [bin_counts bin_values] = imhist(a, num_bins);
 
    % compute standard deviation
    total_pixels = numel(a);
    
    sum_of_pixels = sum(bin_counts .* bin_values);
    mean_pixel = sum_of_pixels / total_pixels;
    
    bin_value_offsets      = bin_values - mean_pixel;
    bin_value_offsets_sqrd = bin_value_offsets .^ 2;
    
    offset_summation = sum( bin_counts .* bin_value_offsets_sqrd);
    s = sqrt(offset_summation / total_pixels);
    
else
    
    % use simple implementation
    if ~isa(a,'double')
        a = double(a);
    end
    s = std(a(:));
    
end