Charger et sauvegarder des images 3D

[Kangourou]

Bonjour,

Cela fait plusieurs fois que je vois passer des questions sur le chargement d'images 3D, ou piles d'images. Il y a quelques temps j'avais bricolé des fonctions pour lire et écrire de telles images. Ces fonctions se basent sur les fonctions imread et imwrite, de la boîte à outils "Image Processing", qu'il vous faudra donc posséder.

Fonction readstack (ou "charge pile") : charge un fichier ou une succession de fichiers. La fonction permet d'importer :

• des séries d'images numérotées en ordre croissant (img00.tif, img01.tif...), en général il suffit de donner la première image pour que toute la pile se charge
• des images tif multi-page
• des images stockées en "raw", dans ce cas il faut spécifier le type de données et la taille du résultat

Code :

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
function img = readstack(fname, varargin)
%READSTACK Read either a list of 2D images (slices), or a 3D image
%
%   Syntax
%   IMG = readstack(FNAME)
%   IMG = readstack(FNAME, INDICES)
%   IMG = readstack(FNAME, TYPE, DIM)
%
%   Description
%
%   IMG = readstack(FNAME)
%   FNAME is the base name of the image.
%   In the case of a bundle image, this is simply the name of the file.
%   In the case of an image stored as several slices, this is the string
%   common to each file with wildcard '#', '?' or '0' at the position of
%   the slices indices.
%   IMG is a 3-dimensional array of type uint8, size: X*Y*Z;
%    or a 4 dimensional array for color images (X*Y*3*Z).
%
%   IMG = readstack(FNAME, INDICES)
%   forces the number of dimensions of resulting array (slices images)
%   FNAME: base filename of images, without end number (string)
%   INDICES: indices of images to put for result. Ex:  [0:39]
%
%   IMG = readstack(FNAME, TYPE, DIM)
%   forces the type and number of dimensions of resulting array.
%   FNAME: name of the single stack file
%   TYPE:  matlab type of data ('uint8', 'double'...)
%   DIM:   size of final stack. Ex : [256 256 50].
%   There is no support for binary raw stacks.
%
%   Examples:
%   IMG = readstack('images???.tif');
%   IMG = readstack('aStack.tif');
%   IMG = readstack('files00.bmp');
%   IMG = readstack('files00.bmp', 5:23);
%   IMG = readstack('files00slices##.tif');
%   IMG = readstack('rawData.raw', 'uint8', [256 256 50]);
%
%   See also:
%   savestack, imread
%
%   ---------
%   author: David Legland, david.legland(at)grignon-dot-inra.fr
%   created the 10/09/2003.

%% import a raw file. 

% First argument is type of image ('uint8' usually), 
% and second argument contains size of image to load.
if length(varargin)>1
    type = varargin{1};
    dim = varargin{2};
    
    img = zeros(dim, type);
    f = fopen(fname, 'r');
    img(:) = fread(f, dim(1)*dim(2)*dim(3), ['*' type]);
    
    return;
end

% check if image stored in one bubndle file or as several stacks
bundle =  sum(ismember('#?', fname))==0;
if exist(fname, 'file')
    bundle = length(imfinfo(fname))>1;
end

if ~bundle    
    % images stored in several 2D files ------------------------
    % -> need to know numbers of slices to read.
        
    %% Compute the base name of the image, by removing '#', '?' or '0'
    
    % characters to replace with numbers
    chars = '#?0';
    
    index = [];    
    for c=1:length(chars)
        % identify number of chars in file name
        n=5;
        while isempty(index) && n>1
            n = n-1;
            index = strfind(fname, repmat(chars(c), [1 n]));
        end

        % escape if special characters has been found
        if  ~isempty(index)
            break;
        end
    end
    
    % In the case of several '00' parts, consider only the last one 
    index = index(length(index));
        
    % create file basename and endname
    len = length(fname);
    basename = fname(1:index-1);
    endname = fname(index+n:len);

    
    %% compute indices of slices to read
    
    if ~isempty(varargin)
        % slice indices are given as parameters
        range=varargin{1};
    else            
        % identify slices to read by detecting last index of slices
        i=0;
        while true
            % check existence of file for given index
            imgname = sprintf(['%s' sprintf('%%0%dd', n) '%s'], basename, i, endname);
            if ~exist(imgname, 'file')
                break;
            end
            i = i+1;
        end
        % read slices from the first one to the last existing one
        range = 0:i-1;
    end
    
    %% Read each slice of the image
    
    % read slice for each range index
    string = [basename sprintf('%%0%dd', n)  endname];
    
    % adapt filename format to windows if needed
    if ispc
        string = strrep(string, '\', '\\');
    end    
    
    % read the first image
    img = imread(sprintf(string, range(1)));
    
    % read each image one after the other
    if length(size(img))==2
        % allocate memory
        img(1, 1, length(range))=0;
        
        % read each gray scale image successively
        for i=2:length(range)
            img(:,:,i) = imread(sprintf(string, range(i)));
        end
    else
        % pre-allocate memory
        img(1, 1, 1, length(range))=0;
        
        % read each color image successively
        for i=2:length(range)
            img(:,:,:,i) = imread(sprintf(string, range(i)));
        end
    end
    
else
    % Image stored in one bundle file --------------------------
    
    % get image dimensions
    info = imfinfo(fname);
    
    % If input argument is found, it is used as the number of slices to
    % read. Anyway, read all the slices.
    range = 1:length(info);
    if ~isempty(varargin)
        range = varargin{1};
    end

    % read each slice of the 3D image
    img = imread(fname, range(1));
    
    if length(size(img))==2         % read gray scale images -----        
        % pre-allocate memory
        img(1, 1, length(range)) = 0;
        
        % add each slice successively
        for i=2:length(range)
            img(:,:,i) = imread(fname, range(i));
        end
    else                            % read color images      -----
        % pre-allocate memory
        img(1, 1, 1, length(range)) = 0;
        
        % add each slice successively
        for i=2:length(range)
            img(:,:,:,i) = imread(fname, range(i));
        end
    end 
end

Et sa petite soeur, la fonction savestack, qui sauve un tableau 3D ou 4D dans une série de fichiers. Pour écrire dans une série de fichier, on passe soit par une chaîne de la forme "img%03d.tif", soit par une succession de caractères "#" ou "?", qui seront remplacés par des chiffres.

Code :

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
function savestack(img, fname, varargin)
%SAVESTACK Save an image stack to a file or a serie of files
%
%   If file name contains '??', '##' or '%0xd' (x being an integer), then
%   the image is saved into a series of files, with increasing index.
%
%   Example :
%   savestack(img, 'imgBundle.tif');        % save as stack
%   savestack(img, 'imgBundle???.tif');     % save as image series
%   savestack(img, 'imgBundle###.tif');     % save as image series
%   savestack(img, 'imgBundle%03d.tif');    % save as image series
%   
%   savestack(..., OPTIONS)
%   use options to write each slice. See imwrite to details.
%
%   See also:
%   readstack, imwrite
%
%   ---------
%   author: David Legland, david.legland(at)grignon.inra.fr
%   created the 10/09/2003.

dim = size(img);

% replace wildcard '??' or '???' by wildcard '%02d' or '%03d'
pos = strfind(fname, '?');
npos = length(pos);
if npos>0
    fname = strrep(fname, repmat('?', [1 npos]), ['%0' num2str(npos) 'd']);
end

% replace wildcard '##' or '###' by wildcard '%02d' or '%03d'
pos = strfind(fname, '#');
npos = length(pos);
if npos>0
    fname = strrep(fname, repmat('#', [1 npos]), ['%0' num2str(npos) 'd']);
end

pos = strfind(fname, '%0');
%if length(strfind(fname, '%'))>0 | slices
if ~isempty(pos)
    % save a series of file, one file per slice of images
    
    % keep only the last '%0xd'
    pos = pos(end);
    
    % extract different parts of the file name
    basename = fname(1:pos-1);
    endname = fname(pos+4:end);
    
    % string format to compute image name
    % -> basename + indSlice + endName
    format = ['%s%0' fname(pos+2) 'd%s'];
    
    % save each slice of the image (gray-scale->3D or color->4D)
    if length(dim)==3
        for i=1:dim(3)
            fileName = sprintf(format, basename, i-1, endname);
            imwrite(img(:,:,i), fileName, ...
                'WriteMode', 'overwrite', varargin{:});
        end
    else
        for i=1:dim(4)
            fileName = sprintf(format, basename, i-1, endname);
            imwrite(img(:,:,:,i), fileName, ...
                'WriteMode', 'overwrite', varargin{:});
        end
    end
    
else
    % save one file containing all slices of image
    if length(dim)==3
        % overwrite existing file
        imwrite(img(:,:,1), fname, varargin{:}, 'WriteMode', 'overwrite');
        
        % append other slices
        for i=2:dim(3)
            imwrite(img(:,:,i), fname, varargin{:}, ...
                'WriteMode', 'append');
        end
    else
        % overwrite existing file
        imwrite(img(:,:,:,1), fname, varargin{:}, 'WriteMode', 'overwrite');

        % append other slices
        for i=1:dim(4)
            imwrite(img(:,:,:,i), fname, varargin{:}, ...
                'WriteMode', 'append');
        end
    end
end

Les deux fonctions permettent de travailler avec des images en niveaux de gris (les images sont des tableaux 3D Ny-Nx-Nz), ou des images couleurs (dans ce cas les images sont des tableaux 4D Ny-Nx-3-Nz).

Bon, je les ai déjà pas mal utilisées, mais un bug est toujours possible. Envoyez les rapports de bug par MP si possible, merci.

A+

[Dut]

Quelques remarques...

Les codes que tu présentes ici ne sont pas identiques à ceux présents dans l'archive zip.

Pour readstack :

Ajoute la ligne suivante au tout début de ton code :

Code :

error(nargchk(1,5,nargin)); % Vérifier le 5

TYPE est déjà un nom de fonction MATLAB, change le ici (et par la suite) :

Code :

type = varargin{1};

Ligne 115, il y a un FOPEN sans FCLOSE correspondant.

Ensuite :

Code :

if length(size(img))==2

s'écrit plus proprement :

Code :

if ndims(img)==2

Le code en commentaire à partir de la ligne 255 ne sert à rien

Code :

string = [basename sprintf('%%0%dd', n)  endname];

vu que tu utilises SPRITNF, s'écrit plus logiquement :

Code :

string = sprintf('%s%%0%dd%s',basename, n, endname);

Code :

1
2
3
	len = length(fname);
    basename = fname(1:index-1);
    endname = fname(index+n:len);

devient simplement :

Code :

1
2
    basename = fname(1:index-1);
    endname = fname(index+n:end);

Code :

index = index(length(index));

devient

Code :

index = index(end);

Pour finir, tu devrais récupérer l'extension des fichiers d'entrée avec FILEPARTS pour la passer ensuite dans IMREAD. Ce sera plus rapide.

Pour savestack

Ajoute la ligne suivante au tout début de ton code :

Code :

error(nargchk(1,4,nargin)); % Vérifier le 4

Il faut aussi ajouter un test sur le format d'image demandé en sortie car tous n'ont pas l'argument 'WriteMode' de IMWRITE (gif, tiff et hdf seulement)

Code :

1
2
3
4
5
pos = strfind(fname, '?');
npos = length(pos);
if npos>0
    fname = strrep(fname, repmat('?', [1 npos]), ['%0' num2str(npos) 'd']);
end

s'écrit plus efficacement :

Code :

1
2
3
4
5
pos = strfind(fname, '?');
if ~isempty(pos)
    npos = length(pos);
    fname = strrep(fname, repmat('?', [1 npos]), ['%0' num2str(npos) 'd']);
end