function proj_geom = gtGeoProjForReconstruction(projvec_sam, omega, vol_center,...
                     bbpos, bboff, detgeo, labgeo, samgeo, recgeo, rectype, vol_center_shift)
% GTGEOPROJFORRECONSTRUCTION Projection coordinates required for
%                            reconstructions.
%
% proj_geom = gtGeoProjForReconstruction(projvec, omega, ...
%             graincenter, bbpos, bboff, labgeo, samgeo, recgeo)
% 
% ----------------------------------------------------------------
%
% Reconstruction geometries
%   'ASTRA_grain' - for a single grain placed in the origin of the 
%                   reconctructed ROI (requires 'graincenter')
%   'ASTRA_full'  - for a reconstruction of the complete sample volume
%                   (does not use 'graincenter')
%   'ASTRA_absorption' - for a reconstruction of the absorption volume
%                   (does not use 'graincenter' and 'projvec')
% 
%
%     INPUT:
%       projvec     = projection vectors (beam direction) in the 
%                     SAMPLE reference (nx3) (Not needed for 'absorption')
%       omega       = omega rotational positions corresponding to 
%                     the projection vectors (in degrees) (nx1)
%       graincenter = grain center in the SAMPLE reference (1x3)
%                     (Only needed for 'grain')
%       bbpos       = bounding box positions on the detector (in pixels)
%                     [corner_U corner_V width height] (nx4)
%       bboff       = additional bounding box U,V offsets relative to
%                     position 'bbpos' (n,2)
%       detgeo      = DETector geometry as in parameters.detgeo
%       labgeo      = LABoratory reference as in parameters.labgeo
%       samgeo      = SAMPLE reference as in parameters.samgeo
%       recgeo      = REConstruction reference as in parameters.recgeo
%       rectype     = <char> reconstruction geometry:
%                     'ASTRA_grain', 'ASTRA_absorption' or 'ASTRA_full'
%       vol_center_shift = shift of grain center in SAMPLE reference (nx3)
%                          accounting for sample movements relative to the
%                          the rotation axis (see also gtMatchGetSampleShifts)
%
%     OUTPUT:
%       proj_geom   = input for the specified reconstruction geometry:
%                       'ASTRA_grain', 'ASTRA_absorption' and 'ASTRA_full':
%                         Parameters defined relative to the sample (nx12):
%                             [projection vector, detector position, ...
%                              detector vector U, detector vector V];
%                       
%

    if isempty(bboff)
        bboff = zeros(size(bbpos, 1), 2);
    end

    if (isempty(detgeo))
        detgeo = gtGeoConvertLegacyLabgeo2Detgeo(labgeo);
    end

    % Bounding box center U,V coordinates
    bbcent_det  = [ ...
        bbpos(:, 1) + (bbpos(:, 3) - 1) / 2 + bboff(:, 1) , ...
        bbpos(:, 2) + (bbpos(:, 4) - 1) / 2 + bboff(:, 2) ];

    if (all([size(omega, 1), size(omega, 2)] == [3 3]))
        % Rotation tensors were given already
        rot_omega = omega;
        ones_omega = ones(size(omega, 3), 1);
    else
        % All functions that use the same omegas, need the same rotation tensors:
        rot_comp_w = gtMathsRotationMatrixComp(labgeo.rotdir', 'col');
        rot_omega = gtMathsRotationTensor(omega, rot_comp_w);
        ones_omega = ones(length(omega), 1);
    end

    if (~exist('vol_center_shift', 'var') || isempty(vol_center_shift))
        vol_center_shift_rec = zeros(numel(ones_omega), 3);
    else
        vol_center_shift_rec = gtGeoSam2Sam(vol_center_shift, samgeo, recgeo, false, false);
    end

    % Bounding box center X,Y,Z coordinates in LAB and REC reference
    bbcent_lab  = gtGeoDet2Lab(bbcent_det, detgeo, 0);
    bbcent_rec  = gtGeoLab2Sam(bbcent_lab, rot_omega, labgeo, recgeo, 0, 0);

    % Detector orientation vectors in REC reference
    detdiru = detgeo.detdiru ./ sqrt(sum(detgeo.detdiru .^ 2)) * detgeo.pixelsizeu;
    detdirv = detgeo.detdirv ./ sqrt(sum(detgeo.detdirv .^ 2)) * detgeo.pixelsizev;

    detdiru_rec = gtGeoLab2Sam(detdiru(ones_omega, :), rot_omega, labgeo, recgeo, true, false);
    detdirv_rec = gtGeoLab2Sam(detdirv(ones_omega, :), rot_omega, labgeo, recgeo, true, false);

    if (strcmpi(rectype, 'ASTRA_absorption'))
        projvec_lab = labgeo.beamdir(ones_omega, :);
        projvec_rec = gtGeoLab2Sam(projvec_lab, rot_omega, labgeo, recgeo, true, false);
    else
        % Normalised projection vector in REC reference
        projvec_rec = gtGeoSam2Sam(projvec_sam, samgeo, recgeo, true, false);
    end
    projvec_rec = gtMathsNormalizeVectorsList(projvec_rec);


    if (strcmpi(rectype, 'ASTRA_grain') || strcmpi(rectype, 'ASTRA_absorption'))
        % Volume center in REC reference
        if (strcmpi(rectype, 'ASTRA_absorption'))
            vol_center_rec = recgeo.orig;
        else
            vol_center_rec = gtGeoSam2Sam(vol_center, samgeo, recgeo, false, false);
        end
        vol_center_rec = vol_center_rec(ones_omega, :) + vol_center_shift_rec;
        bbcent_rec = bbcent_rec - vol_center_rec;
    else
        % case {'ASTRA_full'}
        bbcent_rec = bbcent_rec - vol_center_shift_rec;
    end

    proj_geom = [projvec_rec, bbcent_rec, detdiru_rec, detdirv_rec];
end