% files containing  the images
filenames = {'240px-StLouisArchMultExpEV-4.72.JPG',...
             '240px-StLouisArchMultExpEV-1.82.JPG',...
             '240px-StLouisArchMultExpEV+1.51.JPG',...
             '240px-StLouisArchMultExpEV+4.09.JPG'}; 

%filenames = {'1280px-StLouisArchMultExpEV-4.72.JPG',...
%             '1280px-StLouisArchMultExpEV-1.82.JPG',...
%             '1280px-StLouisArchMultExpEV+1.51.JPG',...
%             '1280px-StLouisArchMultExpEV+4.09.JPG'}; 

imagename = 'hrd_image.png'; % file where the images will be saved

Y = [];
for i = 1:length(filenames),
    % load one of the images
    X = imread(filenames{i});
    
    % convert it in gray levels
    X2 = rgb2gray(X);

    % store it in the multidimensional matrix (in double)
    Y(:,:,i) = X2;
end;

% visualization
for i = 1:size(Y,3),
    figure;
    imshow(uint8(Y(:,:,i)));
end;

% composition by average
Z = mean(Y, 3);
figure;
imshow(uint8(Z));

% save the resulting image
imwrite(uint8(Z), imagename); 