Skip to content
Snippets Groups Projects
gtFwdSimCheckSpotInRawImages.m 2.65 KiB
function [found_intensity, spotInfo] = gtFwdSimCheckSpotInRawImages(uvw, spotsCommProps, parameters)
% CHECKSPOTINRAWIMAGES
%   [found_intensity, spotInfo] = gtFwdSimCheckSpotInRawImages(uvw, spotsCommProps, parameters)
%
% uvw            : spot position at 'u', 'v', 'w'
% spotsCommProps : diffraction spots common properties (Xsize, Ysize,
%                  stackUSize, stackVSize)
% parameters     : parameters.mat

    omega_step = parameters.labgeo.omstep;

    % number of images summed up when searching raw data for a missing spot
    imagerange = round(parameters.fsim.omegarange / omega_step);

    % in case bb_size for search BoundingBox (used for check of
    % completeness) is empty, use the same size as the difspot average size +
    % 5 percent
    if isempty(parameters.fsim.bb_size)
        bbhsize = spotsCommProps.stackUSize;
        bbvsize = spotsCommProps.stackVSize;
    else
        bbhsize = parameters.fsim.bb_size(1);
        bbvsize = parameters.fsim.bb_size(2);
    end

    % Adjust search bbox in order to avoid reading information
    % stemming from regions outside the image border
    bbox_ustart = round(uvw(1) - bbhsize/2);
    bbox_vstart = round(uvw(2) - bbvsize/2);
    bbox_uend   = min(parameters.acq.xdet, bbox_ustart + bbhsize);
    bbox_vend   = min(parameters.acq.ydet, bbox_vstart + bbvsize);
    bb(1) = max(1, bbox_ustart);
    bb(2) = max(1, bbox_vstart);
    bb(3) = (bbox_uend - bbox_ustart); %was +1
    bb(4) = (bbox_vend - bbox_vstart); %was +1
    start_image = max(0, round(uvw(3)) - imagerange);
    end_image   = min(round(uvw(3)) + imagerange, 2*parameters.acq.nproj - 1);

    % read in small image stack at predicted position read ROI from
    % raw images - meanVals are the mean value for each image in the stack
    [raw_stack, meanVals] = gtGetRawRoi(start_image, end_image, parameters.acq, bb);
    [intensity, ~] = max(meanVals);

    % here we force the intensity to be the same in all spots
    volume = spotsCommProps.Xsize * spotsCommProps.Xsize * spotsCommProps.Ysize;
    spot = sum(raw_stack, 3);
    spot = spot / sum(spot(:)) * volume;

    % Checks if found intensity in FULL images at predicted position
    found_intensity = (intensity > parameters.fsim.thr_check);

    % save info of found spot
    spotInfo.raw_stack       = raw_stack;
    spotInfo.meanVals        = meanVals;
    spotInfo.bb              = bb;
    spotInfo.start_image     = start_image;
    spotInfo.end_image       = end_image;
    spotInfo.imagerange      = imagerange;
    spotInfo.found_intensity = found_intensity;
    spotInfo.spot            = spot;
    spotInfo.difblob         = permute(raw_stack,[2 1 3]);
    spotInfo.volume          = volume;
end