%% estimate e, phi, and pr(phi > 0) 
%
% Pools data for 66-68 year olds around 1990 and around age 70, for 90-99
function [ests, fit70, fit90, estBS, bunchBS] = estPooled(bw, opt90, opt70, density, optOuter,p0, name, scale, method, print)
    
    % Load the data around 1990
    yMin = opt90.yMin;
    yMax = opt90.yMax;
    t0 = opt90.t0;
    t1Pre = opt90.t1Pre;
    t1Post = opt90.t1Post;
    ages = opt90.ages;
    [data90, bs90] = setDynamicData(bw, yMin, yMax ,t0, t1Pre, t1Post, ...
        density, ages, scale, method);
    data90.wMat = eye(yMax-yMin+1);
    data90.opt = opt90.opt;
    
    % load the data around age 70
    ages = opt70.ages;
    t0 = opt70.t0;
    t1Pre = opt70.t1Pre;
    t1Post = opt70.t1Post;
    agesD = opt70.agesD;
    [data70, bs70] = setData6970(bw, ages,t0,t1Pre,t1Post,density, agesD, scale);
    data70.wMat = eye(length(ages));
    data70.method = method;
    
    % clean up some to get dynamics properly set
    data70.t0 = ones(size(data70.t1)).*data70.t0;
    data70.zStar = ones(size(data70.t1)).*data70.zStar;
    data70.last = ages == 69;
    data70.ts = ages;
    data70.ts(ages<70) = ages(ages<70)-64;
    data70.ts(ages>=70) = ages(ages>=70)-69;
    data70.opt = opt70.opt;
    
    % put both into a single structure - confusing naming conventions
    data.data70 = data70;
    data.data90 = data90;
    
    % estimation
    [ests, ~, flag] = fminsearch(@(p) calcGmmPooled(p, data), p0 , optOuter);
       
    % Return bootstrap distribution of estimates, if requested
    if (nargout>=4)
        
        nBoot = size(bs70, 2);
        estBS = zeros(nBoot, length(p0));
        dataBS.data70 = data70;
        dataBS.data90 = data90;
        
        for bs = 1:nBoot
           dataBS.data70.bunch = bs70(:, bs);
           dataBS.data90.bunch = bs90(:, bs);
           
           [estsT, ~, flag] = fminsearch(@(p) ...
               calcGmmPooled(p, dataBS), ests, optOuter);
           estBS(bs, :) = estsT;
%            fprintf('.');
%            if bs == 50*round(bs/50)
%                fprintf('\n');
%            end
           if (flag ~= 1)
             fprintf('\n Error: bad flag for %g \n', bs)
           end

        end
        
        bunchBS = [bs70 ; bs90];

    end
    
    % print results 
    mprs = exp(ests(3:end))./(1+exp(ests(3:end)));
    %cprs = cumprod(mprs);
    if print
        fprintf(' %s \\\\ \n', name);
        fprintf(' Elasticity & %4.3f \\\\ \n', ests(1));
        fprintf(' Adjustment cost & %4.1f \\\\ \n', 1e3*ests(2));
        for ii = 1:length(mprs)
            fprintf(' $\\pi_%g$ & %4.3f \\\\ \n', ii, mprs(ii));
        end
    end 
    assert (flag == 1);
    
    % goodness of fit
        % calculate bunching and error for 1990
        % Calculate bounds
        [z1Bar, z0Bar, uz, flag] = calcCutoffs(ests, data90);
        assert (flag == 1);

        % calculate bunching (section 5.2 of dynamic_bunching_model)
        bSim = calcBunching(ests, data90, z1Bar, z0Bar, uz);

        % return objective function
        fit90 = bSim;
        
        if print 
            error90 = bSim - data90.bunch;

            fprintf('\\\\ \\multicolumn{3}{c}{Panel A: Fit for 1990} \\\\ \\midrule \n');
            for ii = 1:length(data90.year);
                fprintf('  %4.0f & %4.3f & %4.3f \\\\ \n', ...
                    data90.year(ii), data90.bunch(ii), bSim(ii));
            end
        end

        % calculate bunching and error for age 70
        % Calculate bounds
        [z1Bar, z0Bar, uz, flag] = calcCutoffs(ests, data70);
        assert (flag == 1);

        % calculate bunching (section 5.2 of dynamic_bunching_model)
        bSim = calcBunching(ests, data70, z1Bar, z0Bar, uz);

        % return objective function
        fit70 = bSim;
        error70 = bSim - data70.bunch;
            
        if print
            fprintf('\\\\ \\multicolumn{3}{c}{Panel B: Fit for age 70} \\\\ \\midrule \n');
            for ii = 1:length(opt70.ages);
                fprintf('  %4.0f & %4.3f & %4.3f \\\\ \n', ...
                    opt70.ages(ii), data70.bunch(ii), fit70(ii));
            end
        end
end
