gtPlaceSubVolume.m

function output = gtPlaceSubVolume(output, input, shift, index, assign_op, use_c_functions)
% GTPLACESUBVOLUME  Analogous to gtPlaceSubImage
%     output = gtPlaceSubVolume(output, input, shift, index, assign_op, use_c_functions)
%     ----------------------------------------------------------------------------------
%     places the input volume in the output vol, with the origin of the input
%     volume at point defined by shift=[x y z] in the output volume.  Any part of the
%     input volume falling outside of the output volume is cropped.  Origin
%     can contain negative coordinates, again anything outside the output
%     volume is cropped.
%     so...  shift=[0 0 0] means that the voxel (a,b,c) in input will be at
%     point (a,b,c) in the output.
%     shift=[2 -2 0] means that voxel (a,b,c) will be at (a+2, b-2, c) in
%     output.
%
%     Note - this used to add on volume to another, rather than place one image
%     in the other.  I have modified so it no longer adds.  Hope this doesn't
%     break anything.
%
%     Note 2 - If index is 0, the volume will be copied as it is, otherwise all
%     the voxels will have value of the index

    if (~exist('use_c_functions', 'var'))
        use_c_functions = true;
    end
%     if (use_c_functions)
%         % Check that the C functions exist!
%         if (~(exist('internal_gtAssignGrainToVol_interf', 'file') == 3 ...
%                 && exist('internal_gtAssignGrainToVol_sum', 'file') == 3 ...
%                 && exist('internal_gtAssignGrainToVol', 'file') == 3 ))
%             use_c_functions = false;
%         end
%     end
    if (~exist('index', 'var'))
        index = 0;
    end
    if (~exist('assign_op', 'var'))
        assign_op = 'assign';
    end

    % Force third dimension to be explicitly given
    [inputSize(1), inputSize(2), inputSize(3)] = size(input);
    [outputSize(1), outputSize(2), outputSize(3)] = size(output);

    % output and input volume limits
    [outLims, inLims] = gtGetVolsIntersectLimits(outputSize, inputSize, shift);

    input = input(  inLims(1, 1):inLims(2, 1), ...
                    inLims(1, 2):inLims(2, 2), ...
                    inLims(1, 3):inLims(2, 3)  );
    lims = [  outLims(1, 1) outLims(2, 1) ...
              outLims(1, 2) outLims(2, 2) ...
              outLims(1, 3) outLims(2, 3)  ];

    if (use_c_functions)
        switch (assign_op)
            case {'zero', 'conflict', 'adaptive'} % Mark overlapping...
                output = internal_gtAssignGrainToVol_interf( output, input, index, lims );
            case 'summed'
                output = internal_gtAssignGrainToVol_sum( output, input, index, lims );
            case {'assign','parent'}
                output = internal_gtAssignGrainToVol( output, input, index, lims );
            otherwise
                error('PLACE:wrong_argument', 'No option for "%s"', assign_op);
        end
    else
        switch (assign_op)
            case {'zero', 'conflict', 'adaptive'}
                % Matlab sugars to do the same that 'internal_gtAssignGrainToVol_interf'
                % does: (please keep them around, in case you break the C function)
                temp_out = output(lims(1):lims(2), lims(3):lims(4), lims(5):lims(6));

                temp_out(temp_out & input) = -1;
                indexes = ((temp_out ==0) & input);
                temp_out(indexes) = input(indexes);

                output(lims(1):lims(2), lims(3):lims(4), lims(5):lims(6)) = temp_out;
            case 'summed'
                % Matlab sugar to do the same that 'internal_gtAssignGrainToVol_sum'
                % does: (please keep them around, in case you break the C function)
                output(lims(1):lims(2), lims(3):lims(4), lims(5):lims(6)) ...
                        = output(lims(1):lims(2), lims(3):lims(4), lims(5):lims(6)) ...
                            + input;
            case {'assign','parent'}
                % Matlab sugar to do the same that 'internal_gtAssignGrainToVol'
                % does: (please keep them around, in case you break the C function)
                output(lims(1):lims(2), lims(3):lims(4), lims(5):lims(6)) = input;

            otherwise
                error('PLACE:wrong_argument', 'No option for "%s"', assign_op);
        end
    end
end