function [ perf ] = frankWolfeOptimizationY( ...
    x, zGt, yGt, yExpected, eventList, longueurs, K, z, y, phi, lambda, mu, nu, params )

%mu is the trade of between the discriminative part and the shape prior
perfZ   = zeros(1, params.niter);
perfY   = zeros(1, params.niter);
jacNoBg = zeros(1, params.niter);
delay   = zeros(1, params.niter);

%If data are in a cell array, we have to stack them
if iscell(x)
    X = cell2mat(x).';
else
    X=x;
end

[N, ~]  = size(X);

% pre-computing heavy stuff
GXTP = computeGXTP(X, lambda);
L = sparse(N, N);

% computing the gradient

params.computeGradFlag = 1;

grad = computeGradient(X, z, N, GXTP);

    if params.silence == 0
        %In that case we pad the Z gradient with 0 otherwise it is not needed
        gradAug = [zeros(1, size(grad, 2)); grad];
    else
        gradAug = grad;
    end
    
    params.grad = grad;

params.computeGradFlag = 0;

obj         = struct('f', [], 'd', [], 't', []);


% computing objective value
f = computeObjective(X, z, y, yExpected, mu, nu, GXTP, params);

% rounding
[z_p, y_p] = rounding(z, K, longueurs, eventList, params);
%     [y_p] = roundingY(z, K, longueurs, eventList, params);

% compare to groundtruth in terms of area under curves (in the Z space)

perfZ0           = lossAssignments(zGt', z_p', K, longueurs, params);
[perfY0, delay0] = lossAlignments(yGt', y_p', K, longueurs, params);

% keeping track of the objective

[res]       = evaluate(z, zGt, z_p, longueurs, K, params);
f10      =  mean(res.f1_per_class);

fprintf('iter=0 f=%-+4.3e ', f);
fprintf('perfZ=%-+4.3e ', perfZ0);
fprintf('perfY=%4.3e', perfY0);
fprintf('meanDelay=%f', delay0);
fprintf('f1=%f', f10);
fprintf('\n ');

tic;

for i = 1:params.niter
    
    % cutting the gradient into clips
    l = mat2cell(gradAug, size(gradAug, 1), longueurs);
    
    %Get the total gradient in Y
    gradY = gradientYFromZ(l, phi, mu, nu, yExpected, y, params);
%      params.computeGradFlag = 1;
% % 
%      dbGradient(z, y, yExpected, gradY, X, mu, nu, phi, GXTP, params) 
%      params.computeGradFlag = 0;

    % performing the Dynamic Programing
    [a, b] = optimizeAY(gradY, phi, params);
    
    if params.silence == 0
        a = a(:, 2:end);
    end
    
    % Linearization duality gap
    
    d = dualityGapY(y, b, gradY);
    
    % getting the optimal step size
    

    params.computeGradFlag = 1;

    gamma = computeGammaOptimalFW(X, z, a, b, GXTP, mu, nu, y, yExpected, params);
%       dbLineSearch(z, y, yExpected, a, b, gamma, X, mu, nu, GXTP, params) 
%     gamma = 2/(i+2);
    % updating z and y
    
    z = (1-gamma) * z + gamma * a;
    y = updateY(gamma, y, b);
    
    

    
    % computing the gradient for the discriminative term
    

    grad = computeGradient(X, z, N, GXTP);
    

    
    if params.silence == 0
        %In that case we pad the Z gradient with 0 otherwise it is not needed
        gradAug = [zeros(1, size(grad, 2)); grad];
    else
        gradAug = grad;
    end
        
     params.computeGradFlag = 0;
%     params.computeGradFlag = 1;
    params.grad = grad;
    
    % computing objective value
    f = computeObjective(X, z, y, yExpected, mu, nu, GXTP, params);
    
    
    % rebuilding the classifiers
    z = getZFromY(y, phi, params);
    w = GXTP * z;
    b = ones(1, N) * (z - X * w) / N;
    
    % rounding
    [z_p, y_p]  = rounding(z, K, longueurs, eventList, params); % z_p is a aggregated matrix , y_p is a cell array
    % [z_p, y_p] = roundingW(x, w, b, yExpected, phi, mu, nu, params);   % rounding with classifier w: un-comment to use this method

    % compare to groundtruth in terms of area under curves (in the Z space)
    
    perfZ(i)               = lossAssignments(zGt', z_p', K, longueurs, params);
    [perfY(i), delay(i)]   = lossAlignments(yGt', y_p', K, longueurs, params);

    % keeping track of the objective
    
    obj(i).f    = f;
    obj(i).d    = d;
    obj(i).t    = toc;
    [res]       = evaluate(z, zGt, z_p, longueurs, K, params);
    f1(i)       =  mean(res.f1_per_class);
    
    % printing the score
    
    fprintf('iter=%3i f=%-+4.3e ', i, f);
    fprintf('dgap=%-+4.3e ', d);
    fprintf('perfZ=%-+4.3e ', perfZ(i));
    fprintf('perfY=%4.3e', perfY(i));
    fprintf('meanDelay=%f', delay(i));
    fprintf('f1=%f', f1(i));
    fprintf('\n ');
    
    
    
end

%Rounding en Y
[y_proundY] = roundingY(y, K, longueurs, eventList, params);
[perfY2, delay2] = lossAlignments(yGt', y_proundY', K, longueurs, params);

yEnd = y;
zEnd = mat2cell(z', K, longueurs);

if params.silence == 0
    silences = 1 - zGt * ones(1, K-1)';
    zGt      = [silences , zGt];
end

zGt = mat2cell(sparse(zGt)', K, longueurs);
params.grad = 0;

%u = randperm(length(y));
%u = u(1:params.numberSaved);

perf.b           = b;
% perf.eventList   = eventList;
perf.w           = w;
%perf.phi         = sparseCell(phi(u));
% perf.zGt         = sparse(zGt);
%perf.yGt         = sparseCell(yGt(u));
%perf.y           = yEnd;
% perf.z           = zEnd;
perf.perfZ       = perfZ;
perf.jac         = jacNoBg;
perf.perfY       = perfY;
perf.delayroundY = delay2;

perf.delay       = delay;
perf.obj         = obj;
perf.lambda      = lambda;
perf.mu          = mu;
params.nu        = nu;
perf.params      = params;

end