source: trunk/src/series/calc_background.m @ 1061

Last change on this file since 1061 was 1061, checked in by g7moreau, 2 years ago
  • Update copyright to year 2019
File size: 9.6 KB
Line 
1%'sub_background': substract background to an image series, used with series.fig
2%------------------------------------------------------------------------
3% function GUI_input=aver_stat(num_i1,num_i2,num_j1,num_j2,Series)
4%
5%OUTPUT
6% GUI_input=list of options in the GUI series.fig needed for the function
7%
8%INPUT:
9%num_i1: series of first indices i (given from the series interface as first_i:incr_i:last_i, mode and list_pair_civ)
10%num_i2: series of second indices i (given from the series interface as first_i:incr_i:last_i, mode and list_pair_civ)
11%num_j1: series of first indices j (given from the series interface as first_j:incr_j:last_j, mode and list_pair_civ )
12%num_j2: series of second indices j (given from the series interface as first_j:incr_j:last_j, mode and list_pair_civ)
13%Series: Matlab structure containing information set by the series interface
14%       .RootPath: path to the image series
15%       .RootFile: root file name
16%       .FileExt: image file extension
17%       .NomType: nomenclature type for file indexing
18%       .NbSlice: %number of slices defined on the interface
19%----------------------------------------------------------------------
20% Method:
21%    calculate the background image by sorting the luminosity of each point
22%     over a sliding sub-sequence of 'nbaver_ima' images.
23%     The luminosity value of rank 'rank' is selected as the
24%     'background'. rank=nbimages/2 gives the median value.  Smaller values are appropriate
25%     for a dense set of particles. The extrem value rank=1 gives the true minimum
26%     luminosity, but it can be polluted by noise.
27% Organization of image indices:
28%     The program is working on a series of images, labelled by two indices i and j, given
29%     by the input matlab vectors num_i1 and num_j1 respectively. In the list, j is the fastest increasing index.
30%     The processing can be done in slices (number nbslice), with bursts of
31%     nbfield2 successive images for a given slice (mode 'multilevel')
32%     In the mode 'volume', nbfield2=1 (1 image at each level)
33
34%=======================================================================
35% Copyright 2008-2019, LEGI UMR 5519 / CNRS UGA G-INP, Grenoble, France
36%   http://www.legi.grenoble-inp.fr
37%   Joel.Sommeria - Joel.Sommeria (A) legi.cnrs.fr
38%
39%     This file is part of the toolbox UVMAT.
40%
41%     UVMAT is free software; you can redistribute it and/or modify
42%     it under the terms of the GNU General Public License as published
43%     by the Free Software Foundation; either version 2 of the license,
44%     or (at your option) any later version.
45%
46%     UVMAT is distributed in the hope that it will be useful,
47%     but WITHOUT ANY WARRANTY; without even the implied warranty of
48%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
49%     GNU General Public License (see LICENSE.txt) for more details.
50%=======================================================================
51
52function GUI_input=calc_background (num_i1,num_i2,num_j1,num_j2,Series)
53
54%------------------------------------------------------------------------
55%requests for the visibility of input windows in the GUI series  (activated directly by the selection in the menu ACTION)
56if ~exist('num_i1','var')
57    GUI_input={'RootPath';'on';...
58        'SubDir';'off';... % subdirectory of derived files (PIV fields), ('on' by default)
59        'RootFile';'on';... %root input file name ('on' by default)
60        'FileExt';'on';... %inputf file extension ('on' by default)
61        'NomType';'on';...%type of file indexing ('on' by default)
62        'NbSlice';'on'; ...%nbre of slices ('off' by default)
63        %'VelTypeMenu';'on';...% menu for selecting the velocity type (civ1,..)('off' by default)
64        %'FieldMenu';'on';...% menu for selecting the velocity field (s) in the input file ('off' by default)
65        %'VelTypeMenu_1';'on';...% menu for selecting the velocity type (civ1,..)('off' by default)
66        %'FieldMenu_1';'on';...% menu for selecting the velocity field (s) in the input file ('off' by default)
67        %'CoordType';...%can use a transform function
68        %'GetObject';...;%can use projection object
69        %'GetMask';...;%can use mask option 
70        %'PARAMETER';'NbSliding';...
71        %'PARAMETER';'VolumeScan';...
72        %'PARAMETER';'RankBrightness';...
73               ''};
74    return %exit the function
75end
76
77%----------------------------------------------------------------
78%% initiate the waitbar
79hseries=guidata(Series.hseries);%handles of the GUI series
80WaitbarPos=get(hseries.waitbar_frame,'Position');
81%-----------------------------------------------------------------
82if iscell(Series.RootPath)
83    msgbox_uvmat('ERROR','This function use only one input image series')
84    return
85end
86
87%% determine input image type
88FileType=[];%default
89MovieObject=[];
90FileExt=Series.FileExt;
91
92if isequal(lower(FileExt),'.avi')
93    hhh=which('mmreader');
94    if ~isequal(hhh,'')&& mmreader.isPlatformSupported()
95        MovieObject=mmreader(fullfile(RootPath,[RootFile FileExt]));
96        FileType='movie';
97    else
98        FileType='avi';
99    end
100elseif isequal(lower(FileExt),'.vol')
101    FileType='vol';
102else
103   form=imformats(FileExt(2:end));
104   if ~isempty(form)% if the extension corresponds to an image format recognized by Matlab
105       if isequal(Series.NomType,'*');
106           FileType='multimage';
107       else
108           FileType='image';
109       end
110   end
111end
112if isempty(FileType)
113    msgbox_uvmat('ERROR',['invalid file extension ' FileExt ': this function only accepts image or movie input'])
114    return
115end
116
117nbslice_i=Series.NbSlice; %number of slices
118siz=size(num_i1);
119nbaver_init=23;%approximate number of images used for the sliding background: to be adjusted later to include an integer number of bursts
120
121filebase=fullfile(Series.RootPath,Series.RootFile);
122dir_images=Series.RootPath;
123nom_type=Series.NomType;
124
125%% create dir of the new image
126% [dir_images,namebase]=fileparts(filebase);
127term='_background';
128
129[pp,subdir_ima]=fileparts(Series.RootPath);
130try
131    mkdir([dir_images term]);
132catch ME
133            msgbox_uvmat('ERROR',ME.message);
134            return
135end
136[xx,msg2] = fileattrib([dir_images term],'+w','g'); %yield writing access (+w) to user group (g)
137if ~strcmp(msg2,'')
138    msgbox_uvmat('ERROR',['pb of permission for ' subdir_ima term ': ' msg2])%error message for directory creation
139    return
140end
141filebase_b=fullfile([dir_images term],Series.RootFile);
142
143%% set processing parameters
144prompt = {'The number of positions (laser slices)';'volume scan mode (Yes/No)';...
145        'the luminosity rank chosen to define the background (0.1=for dense particle seeding, 0.5 (median) for sparse particles'};
146dlg_title = ['get (slice by slice) a  background , result in subdir ' subdir_ima term];
147num_lines= 3;
148def     = {'1';'No';'0.1'};
149answer = inputdlg(prompt,dlg_title,num_lines,def);
150vol_test=answer{2};
151if isequal(vol_test,'Yes')
152    nbfield2=1;%case of volume: no consecutive series at a given level
153    nbslice_i=siz(1);%number of slices
154else
155    nbfield2=siz(1); %nb of consecutive images at each level(burst)
156%     if siz(2)>1
157%        nbslice_i=str2num(answer{1})/(num_i1(1,2)-num_i1(1,1));% number of slices
158%     else
159        nbslice_i=1;
160%     end
161%     if ~isequal(floor(nbslice_i),nbslice_i)
162%         msgbox_uvmat('ERROR','the number of slices must be a multiple of the i increment')
163%         return
164%     end
165end
166lengthtot=siz(1)*siz(2);
167nbfield=floor(lengthtot/(nbfield2*nbslice_i));%total number of i indexes (adjusted to an integer number of slices)
168nbfield_slice=nbfield*nbfield2;% number of fields per slic
169rank=floor(str2num(answer{3})*nbfield_slice);
170if rank==0
171    rank=1;%rank selected in the sorted image series
172end
173
174%% prealocate memory for the background
175first_image=name_generator(filebase,num_i1(1),num_j1(1),Series.FileExt,Series.NomType);
176Afirst=read_image(first_image,FileType,num_i1(1),MovieObject);
177[npy,npx]=size(Afirst);
178try
179    Ak=zeros(npy,npx,nbfield_slice,'uint16'); %prealocate memory
180catch ME
181    msgbox_uvmat('ERROR',ME.message)
182    return
183end
184
185
186%MAIN LOOP ON SLICES
187%
188% nbfield_slice
189% nbslice_i
190% for islice=1:nbslice_i
191%     %% select the series of image indices at the level islice
192%     for ifield=1:nbfield
193%         for iburst=1:nbfield2
194%             indselect(iburst,ifield)=((ifield-1)*nbslice_i+(islice-1))*nbfield2+iburst;
195%         end
196%     end
197%     
198    %% read the first series of nbaver_ima images and sort by luminosity at each pixel
199    for ifile = 1:nbfield_slice
200%         ifile=indselect(ifield);
201        filename=name_generator(filebase,num_i1(ifile),num_j1(ifile),Series.FileExt,Series.NomType)
202        try
203            Aread=read_image(filename,FileType,num_i1(ifile),MovieObject);
204        catch ME
205            msgbox_uvmat('ERROR',ME.message)
206            return
207        end
208        Ak(:,:,ifile)=Aread;
209            %finish the waitbar
210        update_waitbar(hseries.waitbar,WaitbarPos,1)
211    end
212    Ak=sort(Ak,3);
213    C=Ak(:,:,rank);
214    [newname]=...
215        name_generator(filebase_b,num_i1(1),num_j1(1),'.png','none') % makes the new file name 
216    imwrite(C,newname,'BitDepth',16); % save the new image
217% end   
218   
219
220
221
222%------------------------------------------------------------------------
223%--read images and convert them to the uint16 format used for PIV
224function A=read_image(filename,type_ima,num,MovieObject)
225%------------------------------------------------------------------------
226%num is the view number needed for an avi movie
227switch type_ima
228    case 'movie'
229        A=read(MovieObject,num);
230    case 'avi'
231        mov=aviread(filename,num);
232        A=frame2im(mov(1));
233    case 'multimage'
234        A=imread(filename,num);
235    case 'image'   
236        A=imread(filename);
237end
238siz=size(A);
239if length(siz)==3;%color images
240    A=sum(double(A),3);
241end
242   
Note: See TracBrowser for help on using the repository browser.