annotate toolboxes/AudioInpaintingToolbox/Solvers/inpaintSignal_IndependentProcessingOfFrames.m @ 178:4ea4badb2266 danieleb

added ramirez dl (to be completed) and MOCOD dictionary update
author Daniele Barchiesi <daniele.barchiesi@eecs.qmul.ac.uk>
date Thu, 17 Nov 2011 11:22:17 +0000
parents 56d719a5fd31
children
rev   line source
ivan@138 1 function [ReconstSignal1 ReconstSignal2] = inpaintSignal_IndependentProcessingOfFrames(problemData,param)
ivan@138 2 %
ivan@138 3 %
ivan@138 4 % Usage:
ivan@138 5 %
ivan@138 6 %
ivan@138 7 % Inputs:
ivan@138 8 % -
ivan@138 9 % -
ivan@138 10 % -
ivan@138 11 % -
ivan@138 12 % -
ivan@138 13 % -
ivan@138 14 % -
ivan@138 15 % -
ivan@138 16 %
ivan@138 17 % Outputs:
ivan@138 18 % -
ivan@138 19 % -
ivan@138 20 % -
ivan@138 21 % -
ivan@138 22 %
ivan@138 23 % Note that the CVX library is needed.
ivan@138 24 %
ivan@138 25 % -------------------
ivan@138 26 %
ivan@138 27 % Audio Inpainting toolbox
ivan@138 28 % Date: June 28, 2011
ivan@138 29 % By Valentin Emiya, Amir Adler, Maria Jafari
ivan@138 30 % This code is distributed under the terms of the GNU Public License version 3 (http://www.gnu.org/licenses/gpl.txt).
ivan@138 31 % ========================================================
ivan@138 32 % Perform Audio De-clipping with overlapping blocks
ivan@138 33 % Synthesis Approach, union of overcomplete DCT dictionary
ivan@138 34 % Date: 14 Apr. 2010
ivan@138 35 % Inputs:
ivan@138 36 % - x: Clipped signal
ivan@138 37 % - ClipMask: Indices of clipped samples
ivan@138 38 % - Optional parameters [and default values]:
ivan@138 39 % - param.N: frame length [256]
ivan@138 40 % - param.frameOverlapFactor: overlap factor between frames [2]
ivan@138 41 % - param.wa: weighting analysis window [@sineWin]
ivan@138 42 % - param.ws: weighting synthesis window [@sineWin]
ivan@138 43 % - param.OMPerr: error threshold to stop OMP iterations [0.001]
ivan@138 44 % - param.sparsityDegree: max number of non-zero components to
ivan@138 45 % stop OMP iterations [param.N/4];
ivan@138 46 % - other fields: see the documentation of UDCT_Dictionary
ivan@138 47 %
ivan@138 48 % Outputs:
ivan@138 49 % ReconstSignal1 - reconstructed signal (all samples generated
ivan@138 50 % from the synthesis model)
ivan@138 51 % ReconstSignal2 - reconstructed signal (only clipped samples are generated
ivan@138 52 % from the synthesis model)
ivan@138 53 %
ivan@138 54 % By Valentin Emiya - SMALL Project, 2010
ivan@138 55 %
ivan@138 56 % ========================================================
ivan@138 57
ivan@138 58 % Check parameters
ivan@138 59 defaultParam.N = 256;
ivan@138 60 defaultParam.OLA_frameOverlapFactor = 4;
ivan@138 61 defaultParam.wa = @wSine;
ivan@138 62 defaultParam.OLA_ws = @wSine;
ivan@138 63 defaultParam.OLA_par_waitingTime_mainProcess = 0.2;
ivan@138 64 defaultParam.OLA_par_waitingTime_thread = 0.2;
ivan@138 65 defaultParam.OLA_par_frameBlockSize = 1;
ivan@138 66 defaultParam.TCPIP_port = 3000;
ivan@138 67 defaultParam.COM_DISP = false;
ivan@138 68 defaultParam.STATE_DISP = false;
ivan@138 69
ivan@138 70 if nargin<2
ivan@138 71 param = defaultParam;
ivan@138 72 else
ivan@138 73 names = fieldnames(defaultParam);
ivan@138 74 for k=1:length(names)
ivan@138 75 if ~isfield(param,names{k}) || isempty(param.(names{k}))
ivan@138 76 param.(names{k}) = defaultParam.(names{k});
ivan@138 77 end
ivan@138 78 end
ivan@138 79 end
ivan@138 80
ivan@138 81 x = problemData.x;
ivan@138 82 ClipMask = find(problemData.IMiss);
ivan@138 83
ivan@138 84 % According to this flag, switch between a parallel multithread processing
ivan@138 85 % and a singlethread processing. This can fasten the computations but does
ivan@138 86 % not affect the results.
ivan@138 87 if param.MULTITHREAD_FRAME_PROCESSING
ivan@138 88 [ReconstSignal1 ReconstSignal2] = multithreadProcessing(x,ClipMask,param);
ivan@138 89 else
ivan@138 90 [ReconstSignal1 ReconstSignal2] = singlethreadProcessing(x,ClipMask,param);
ivan@138 91 end
ivan@138 92
ivan@138 93 return;
ivan@138 94
ivan@138 95 function [ReconstSignal1 ReconstSignal2] = singlethreadProcessing(x,ClipMask,param)
ivan@138 96 % ========================================================
ivan@138 97 % Overlap-and-add processing of a signal with missing samples.
ivan@138 98 % Decomposition into overlapping frames, processing of each
ivan@138 99 % frame independently and OLA reconstruction.
ivan@138 100 % Date: 01 Jun. 2010
ivan@138 101 % Inputs:
ivan@138 102 % - x: Clipped signal
ivan@138 103 % - ClipMask: Indices of clipped samples
ivan@138 104 % - Optional parameters [and default values]:
ivan@138 105 % - param.N: frame length [256]
ivan@138 106 % - param.inpaintFrame: function handle for inpainting a frame
ivan@138 107 % [@inpaintFrame_OMP]
ivan@138 108 % - param.OLA_frameOverlapFactor: overlap factor between frames [2]
ivan@138 109 % - param.wa: weighting analysis window [@sineWin]
ivan@138 110 % - param.OLA_ws: weighting synthesis window [@sineWin]
ivan@138 111 % - param.SKIP_CLEAN_FRAMES: flag to skip frames with no
ivan@138 112 % missing samples [true]
ivan@138 113 % - other fields: see the documentation the inpainting method
ivan@138 114 %
ivan@138 115 % Outputs:
ivan@138 116 % ReconstSignal1 - reconstructed signal (all samples generated
ivan@138 117 % from the synthesis model)
ivan@138 118 % ReconstSignal2 - reconstructed signal (only clipped samples are generated
ivan@138 119 % from the synthesis model)
ivan@138 120 %
ivan@138 121 % By Amir Adler, Maria Jafari, Valentin Emiya - SMALL Project, 2010
ivan@138 122 %
ivan@138 123 % ========================================================
ivan@138 124
ivan@138 125 % Check parameters
ivan@138 126 defaultParam.N = 256;
ivan@138 127 defaultParam.inpaintFrame = @inpaintFrame_OMP;
ivan@138 128 defaultParam.OLA_frameOverlapFactor = 2;
ivan@138 129 defaultParam.wa = @wSine;
ivan@138 130 defaultParam.ws = @wSine;
ivan@138 131 defaultParam.SKIP_CLEAN_FRAMES = true;
ivan@138 132
ivan@138 133 if nargin<3
ivan@138 134 param = defaultParam;
ivan@138 135 else
ivan@138 136 names = fieldnames(defaultParam);
ivan@138 137 for k=1:length(names)
ivan@138 138 if ~isfield(param,names{k}) || isempty(param.(names{k}))
ivan@138 139 param.(names{k}) = defaultParam.(names{k});
ivan@138 140 end
ivan@138 141 end
ivan@138 142 end
ivan@138 143 % if ~isfield(param,'D')
ivan@138 144 % param.D = param.D_fun(param);
ivan@138 145 % end
ivan@138 146
ivan@138 147 bb=param.N; % block size
ivan@138 148
ivan@138 149 % modify signal length to accommodate integer number of blocks
ivan@138 150 L=floor(length(x)/bb)*bb;
ivan@138 151 x=x(1:L);
ivan@138 152 ClipMask(ClipMask>L) = [];
ivan@138 153
ivan@138 154 % Extracting the signal blocks with 50% overlap
ivan@138 155 Ibegin = (1:bb/param.OLA_frameOverlapFactor:length(x)-bb);
ivan@138 156 if Ibegin(end)~=L-bb+1
ivan@138 157 Ibegin(end+1) = L-bb+1;
ivan@138 158 end
ivan@138 159 Iblk = ones(bb,1)*Ibegin+(0:bb-1).'*ones(size(Ibegin));
ivan@138 160 wa = param.wa(bb); % analysis window
ivan@138 161 xFrames=diag(wa)*x(Iblk);
ivan@138 162
ivan@138 163 % Generating the block mask
ivan@138 164 Mask=ones(size(x));
ivan@138 165 Mask(ClipMask)=0;
ivan@138 166 blkMask=Mask(Iblk);
ivan@138 167
ivan@138 168 % Declipping the Patches
ivan@138 169 [n,P]=size(xFrames);
ivan@138 170
ivan@138 171 if ~isdeployed
ivan@138 172 h=waitbar(0,'Processing each frame...');
ivan@138 173 end
ivan@138 174 Reconst = zeros(n,P);
ivan@138 175 co = zeros(512,P);
ivan@138 176 for k=1:1:P,
ivan@138 177 if ~isdeployed
ivan@138 178 waitbar(k/P);
ivan@138 179 end
ivan@138 180 if param.SKIP_CLEAN_FRAMES && all(blkMask(:,k))
ivan@138 181 continue
ivan@138 182 end
ivan@138 183 frameProblemData.x = xFrames(:,k);
ivan@138 184 frameProblemData.IMiss = ~blkMask(:,k);
ivan@138 185
ivan@138 186 [Reconst(:,k)]= ...
ivan@138 187 param.inpaintFrame(frameProblemData,param);
ivan@138 188
ivan@138 189 end;
ivan@138 190 if ~isdeployed
ivan@138 191 close(h);
ivan@138 192 end
ivan@138 193
ivan@138 194 % Overlap and add
ivan@138 195
ivan@138 196 % The completly reconstructed signal
ivan@138 197 ReconstSignal1 = zeros(size(x));
ivan@138 198 ws = param.OLA_ws(bb); % synthesis window
ivan@138 199 wNorm = zeros(size(ReconstSignal1));
ivan@138 200 for k=1:size(Iblk,2)
ivan@138 201 ReconstSignal1(Iblk(:,k)) = ReconstSignal1(Iblk(:,k)) + Reconst(:,k).*ws(:);
ivan@138 202 wNorm(Iblk(:,k)) = wNorm(Iblk(:,k)) + ws(:).*wa(:);
ivan@138 203 end
ivan@138 204 ReconstSignal1 = ReconstSignal1./wNorm;
ivan@138 205
ivan@138 206 % Only replace the clipped samples with the reconstructed ones
ivan@138 207 ReconstSignal2=x;
ivan@138 208 ReconstSignal2(ClipMask)=ReconstSignal1(ClipMask);
ivan@138 209
ivan@138 210 return;
ivan@138 211
ivan@138 212 function [ReconstSignal1 ReconstSignal2] = multithreadProcessing(x,ClipMask,param)
ivan@138 213 % Send parameters to the threads
ivan@138 214 % initParamFilename = [param.OLA_par_threadDir 'par_param.mat'];
ivan@138 215 % save(initParamFilename,'param');
ivan@138 216
ivan@138 217 bb=param.N; % block size
ivan@138 218
ivan@138 219 % modify signal length to accommodate integer number of blocks
ivan@138 220 L=floor(length(x)/bb)*bb;
ivan@138 221 x=x(1:L);
ivan@138 222 ClipMask(ClipMask>L) = [];
ivan@138 223
ivan@138 224 % Extracting the signal blocks with 50% overlap
ivan@138 225 Ibegin = (1:round(bb/param.OLA_frameOverlapFactor):length(x)-bb);
ivan@138 226 if Ibegin(end)~=L-bb+1
ivan@138 227 Ibegin(end+1) = L-bb+1;
ivan@138 228 end
ivan@138 229 Iblk = ones(bb,1)*Ibegin+(0:bb-1).'*ones(size(Ibegin));
ivan@138 230 wa = param.wa(bb); % analysis window
ivan@138 231 xFrames=diag(wa)*x(Iblk);
ivan@138 232
ivan@138 233 % Generating the block mask
ivan@138 234 Mask=ones(size(x));
ivan@138 235 Mask(ClipMask)=0;
ivan@138 236 blkMask=Mask(Iblk);
ivan@138 237
ivan@138 238 % Declipping the Patches
ivan@138 239 if ~isdeployed && false
ivan@138 240 h=waitbar(0,'Processing each frame...');
ivan@138 241 end
ivan@138 242 [n,P]=size(xFrames);
ivan@138 243 Reconst = NaN(n,P);
ivan@138 244 % initializedThreads = [];
ivan@138 245
ivan@138 246 % find the block of frames to process
ivan@138 247 k_lists = {};
ivan@138 248 kTrame = 1;
ivan@138 249 while kTrame<=P
ivan@138 250 k_list = zeros(param.OLA_par_frameBlockSize,1);
ivan@138 251 ind = 0;
ivan@138 252 while ind<param.OLA_par_frameBlockSize && kTrame<=P
ivan@138 253 if param.SKIP_CLEAN_FRAMES && all(blkMask(:,kTrame))
ivan@138 254 kTrame=kTrame+1;
ivan@138 255 continue
ivan@138 256 end
ivan@138 257 ind = ind+1;
ivan@138 258 k_list(ind) = kTrame;
ivan@138 259 kTrame=kTrame+1;
ivan@138 260 end
ivan@138 261 if ind==0
ivan@138 262 break;
ivan@138 263 end
ivan@138 264 k_lists{end+1} = k_list(1:ind);
ivan@138 265 end
ivan@138 266 k_list_all = cell2mat(k_lists');
ivan@138 267
ivan@138 268 % Create a server
ivan@138 269 serverSocket = createServer(param);
ivan@138 270 % Definition of the client states
ivan@138 271 stateDef;
ivan@138 272
ivan@138 273 param.COM_DISP = false;
ivan@138 274
ivan@138 275 kBlock=1;
ivan@138 276 initializedClientIDs = [];
ivan@138 277 while any(isnan(Reconst(1,k_list_all)))
ivan@138 278 if ~isdeployed && false
ivan@138 279 waitbar(sum(~isnan(Reconst(1,k_list_all)))/length(k_list_all));
ivan@138 280 end
ivan@138 281
ivan@138 282 % Wait for a new client
ivan@138 283 currentClient = waitClient(serverSocket);
ivan@138 284
ivan@138 285 % Receive client state
ivan@138 286 clientIDState = readData(currentClient,param.COM_DISP);
ivan@138 287 clientID = clientIDState(1);
ivan@138 288 clientState = clientIDState(2);
ivan@138 289 switch clientState
ivan@138 290 case INIT
ivan@138 291 if param.STATE_DISP
ivan@138 292 fprintf('INIT %d\n',clientID);
ivan@138 293 end
ivan@138 294 if 0
ivan@138 295 sendData(currentClient,initParamFilename,param.COM_DISP);
ivan@138 296 else
ivan@138 297 sendData(currentClient,param,param.COM_DISP);
ivan@138 298 end
ivan@138 299 initializedClientIDs(end+1) = clientID;
ivan@138 300 case FREE
ivan@138 301 if ~ismember(clientID,initializedClientIDs)
ivan@138 302 sendData(currentClient,INIT_ORDER,param.COM_DISP); % INIT
ivan@138 303 elseif kBlock<=length(k_lists)
ivan@138 304 k_list = k_lists{kBlock};
ivan@138 305 if param.STATE_DISP
ivan@138 306 fprintf('TO PROCESS %d:',clientID);
ivan@138 307 arrayfun(@(x)fprintf(' %d',x),k_list);
ivan@138 308 fprintf('\n');
ivan@138 309 end
ivan@138 310 sendData(currentClient,k_list,param.COM_DISP);
ivan@138 311 sendData(currentClient,xFrames(:,k_list),param.COM_DISP);
ivan@138 312 sendData(currentClient,find(blkMask(:,k_list)),param.COM_DISP);
ivan@138 313 kBlock = kBlock+1;
ivan@138 314 else
ivan@138 315 if param.STATE_DISP
ivan@138 316 fprintf('WAIT %d\n',clientID);
ivan@138 317 end
ivan@138 318 sendData(currentClient,WAIT_ORDER,param.COM_DISP); % no data to process
ivan@138 319 end
ivan@138 320 case PROCESSED
ivan@138 321 processed_k_list = readData(currentClient,param.COM_DISP);
ivan@138 322 y = readData(currentClient,param.COM_DISP);
ivan@138 323 y = reshape(y,[],length(processed_k_list));
ivan@138 324 if param.STATE_DISP
ivan@138 325 fprintf('PROCESSED %d:',clientID);
ivan@138 326 arrayfun(@(x)fprintf(' %d',x),processed_k_list);
ivan@138 327 fprintf('\n');
ivan@138 328 end
ivan@138 329 if ~isempty(processed_k_list)
ivan@138 330 Reconst(:,processed_k_list) = y;
ivan@138 331 end
ivan@138 332 otherwise
ivan@138 333 error('switch:UndefinedClientState','Undefined client state');
ivan@138 334 end
ivan@138 335
ivan@138 336 closeClientConnection(currentClient);
ivan@138 337 end;
ivan@138 338
ivan@138 339 % Close the server
ivan@138 340 closeServer(serverSocket);
ivan@138 341
ivan@138 342 if ~isdeployed && false
ivan@138 343 close(h);
ivan@138 344 end
ivan@138 345
ivan@138 346 % Overlap and add
ivan@138 347
ivan@138 348 % The completly reconstructed signal
ivan@138 349 ReconstSignal1 = zeros(size(x));
ivan@138 350 ws = param.OLA_ws(bb); % synthesis window
ivan@138 351 wNorm = zeros(size(ReconstSignal1));
ivan@138 352 for k=1:size(Iblk,2)
ivan@138 353 ReconstSignal1(Iblk(:,k)) = ReconstSignal1(Iblk(:,k)) + Reconst(:,k).*ws(:);
ivan@138 354 wNorm(Iblk(:,k)) = wNorm(Iblk(:,k)) + ws(:).*wa(:);
ivan@138 355 end
ivan@138 356 ReconstSignal1 = ReconstSignal1./wNorm;
ivan@138 357
ivan@138 358 % Only replace the clipped samples with the reconstructed ones
ivan@138 359 ReconstSignal2=x;
ivan@138 360 ReconstSignal2(ClipMask)=ReconstSignal1(ClipMask);
ivan@138 361
ivan@138 362 killClients(param);
ivan@138 363
ivan@138 364 return;
ivan@138 365
ivan@138 366 % ========================================================
ivan@138 367 % ========================================================
ivan@138 368
ivan@138 369
ivan@138 370
ivan@138 371 % ========================================================
ivan@138 372 % ========================================================