gtReconstructTestCase.m

function [result_viewer, result, post_result, algo] = gtReconstructTestCase(test_data, det_index, do_plot)

    if (~exist('det_index', 'var') || isempty(det_index))
        det_index = 1;
    end
    if (~exist('do_plot', 'var') || isempty(do_plot))
        do_plot = false;
    end

    rec_opts = gtReconstruct6DGetParamenters(test_data.parameters);

    if (~isfield(test_data, 'gv'))
        test_data.gv.dm = gtDefDmvol2Gvdm(test_data.dmvol);
    end
    gvdm = test_data.gv.dm(:, test_data.gv.used_ind);
%     gvcs = test_data.gv.cs(:, test_data.gv.used_ind);
%     if (~isfield(test_data, 'intvol'))
%         size_int_vol = [ ...
%             size(test_data.dmvol, 1), ...
%             size(test_data.dmvol, 2), ...
%             size(test_data.dmvol, 3) ];
%         test_data.intvol = ones(size_int_vol, 'like', test_data.dmvol);
%     end
%     gvpow = test_data.intvol(test_data.gv.used_ind);

    cryst = test_data.parameters.cryst(test_data.grain.phaseid);
    symm = gtCrystGetSymmetryOperators(cryst.crystal_system, cryst.spacegroup);
    gtGetMaxDisorientation(gvdm, symm, 'diameter_average');
    gtGetMaxDisorientation(gvdm, symm, 'diameter_zero');

    if (do_plot)
        f = figure();
        ax = axes('parent', f);
        hold(ax, 'on')
        scatter3(ax, gvdm(1, 1:10:end)', gvdm(2, 1:10:end)', gvdm(3, 1:10:end)', 20);
    end

    sampler = GtOrientationSampling(test_data.parameters, test_data.grain, 'detector_index', det_index(1));

    if (do_plot)
        % Testing ability of GtOrientationSampling to detect the correct
        % orientation-space boundingbox
        sampler.make_grid_resolution_estim_ODF('cubic', ...
            rec_opts.ospace_resolution, 1, rec_opts.max_grid_edge_points);

        scatter3(ax, sampler.R_vectors{1}(:, 1), sampler.R_vectors{1}(:, 2), sampler.R_vectors{1}(:, 3), 20, 'r');
        drawnow();
        hold(ax, 'on')
    end

    % Using the optimal orientation-space boundingbox
    sampler.make_grid_resolution('cubic', rec_opts.ospace_resolution, ...
        gvdm, rec_opts.ospace_oversize, rec_opts.max_grid_edge_points);

    if (do_plot)
%         ospace_bb = [min(gvdm, [], 2) - 0.002, max(gvdm, [], 2) + 0.002];
%         sampler.make_grid('cubic', 3, ospace_bb, rec_opts.ospace_oversize);

        scatter3(ax, sampler.R_vectors{1}(:, 1), sampler.R_vectors{1}(:, 2), sampler.R_vectors{1}(:, 3), 20, 'g');
        hold(ax, 'off')
        drawnow();
    end

    if (rec_opts.ospace_super_sampling > 1)
        sampler.make_supersampling_simple_grid([1 2 3], rec_opts.ospace_super_sampling);
    end

%     ref_gr = sampler.get_reference_grain();
%     proj = ref_gr.proj(det_index);
% 
%     spacing = mean([proj.vol_size_y, proj.vol_size_x, proj.vol_size_z]) * (rec_opts.rspace_oversize - 1);
%     volume_size = ceil([proj.vol_size_y, proj.vol_size_x, proj.vol_size_z] + spacing);
% 
%     if (rec_opts.volume_downscaling > 1)
%         volume_size = ceil(volume_size / rec_opts.volume_downscaling);
%     end
% 
%     r_voxsize = test_data.parameters.recgeo(det_index).voxsize * rec_opts.volume_downscaling;
%     half_volume_size = (volume_size - 1) / 2 .* r_voxsize;
%     rspace_grid_points = [ref_gr.center - half_volume_size, ref_gr.center + half_volume_size, r_voxsize];
% 
%     odf_6D = gtGetODF6DFromGvdm(gvdm, gvcs, sampler.lattice.gr, rspace_grid_points, gvpow, 'linear');
% 
%     vols_6D = reshape(permute(odf_6D, [4 5 6 1 2 3]), volume_size(1), volume_size(2), volume_size(3), []);
%     vols_6D = single(vols_6D);
%     vols = cell(size(vols_6D, 4), 1);
%     for ii_v = 1:size(vols_6D, 4)
%         vols{ii_v} = vols_6D(:, :, :, ii_v);
%     end

    algo = gtReconstruct6DLaunchAlgorithm(sampler, rec_opts, ...
        test_data.parameters, 'det_index', det_index);

    vols = algo.getCurrentSolution();

    or_sizes = sampler.get_orientation_sampling_size();

    [avg_R_vecs, avg_R_vecs_int, stddev_R_vecs] = sampler.getAverageOrientations(vols);
    avg_R_vec = sampler.getAverageOrientation(vols);
    s_g_odf = reshape(sampler.getODF(vols), or_sizes);

    vol_size = size(avg_R_vecs_int);
    shift = gtFwdSimComputeVolumeShifts(test_data.grain.proj, test_data.parameters, vol_size, det_index(1));

    [kam, gam] = gtDefComputeKernelAverageMisorientation(avg_R_vecs, avg_R_vecs_int);
    [igm, gos] = gtDefComputeIntraGranularMisorientation(avg_R_vecs, avg_R_vecs_int, 'R_vector', test_data.grain.R_vector);

    % Restoring initial volume size (depending on the rounding)
    if (rec_opts.volume_downscaling > 1)
        fprintf('Expanding volumes..')
        c = tic();
        avg_R_vecs_int = gtMathsUpsampleVolume(avg_R_vecs_int, rec_opts.volume_downscaling);
        avg_R_vecs = gtMathsUpsampleVolume(avg_R_vecs, rec_opts.volume_downscaling);
        stddev_R_vecs = gtMathsUpsampleVolume(stddev_R_vecs, rec_opts.volume_downscaling);
        kam = gtMathsUpsampleVolume(kam, rec_opts.volume_downscaling);
        igm = gtMathsUpsampleVolume(igm, rec_opts.volume_downscaling);

        vols = gtMathsUpsampleVolume(vols, rec_opts.volume_downscaling);
        fprintf('\b\b: Done (%f seconds).\n', toc(c))
    end

    result = gt6DAssembleResult(test_data, sampler, vols, det_index(1));

    result.ODF6D = struct( ...
        'voxels_avg_R_vectors', {avg_R_vecs}, ...
        'intensity', {avg_R_vecs_int}, ...
        'shift', {shift}, ...
        'R_vectors', {sampler.get_R_vectors()}, ...
        'voxels_stddev_R_vectors', {stddev_R_vecs}, ...
        'single_grain_ODF', {s_g_odf}, ...
        'single_grain_avg_R_vector', {avg_R_vec}, ...
        'kernel_average_misorientation', {kam}, ...
        'intra_granular_misorientation', {igm} );

    post_result = gt6DPostProcessOrientationSpread(test_data, result, true);

    [proj_blobs, proj_spots] = algo.getProjectionOfCurrentSolution();

    post_result.ODF6D = struct( ...
        'orientation_volumes', {vols}, ...
        'fwd_projected_blobs', {proj_blobs}, ...
        'fwd_projected_spots', {proj_spots}, ...
        'grain_average_misorientation', {gam}, ...
        'grain_orientation_spread', {gos}, ...
        'compute_statistics', {algo.get_statistics()} );

    result_viewer = GtOrientationResultView(test_data, result, post_result, ...
        'f_title', sprintf('Result Browser: algo "%s", det_id:%s"', ...
        test_data.parameters.rec.grains.algorithm, sprintf(' %d', det_index)));
end