function fixedStatus=fixFlatMesh1File(inFileName,outFileName,grayFile)
% function fixedStatus=fixFlatMesh1File(inFileName,outFileName)
% mrFlatMesh v1.0 had a nasty bug : the 3dGlocs of layer 1 gnodes were
% wrong by -1 voxel in all directions (i.e. a node at [1,2,3] was 
% given a 3D location of [0,1,2]
% This has been fixed in v1.5 upwards but old flat.mat files need to be fixed.
% This routine reads in a flat.mat file and tries to fix it by adding 1 to all the layer 1 gLocs3ds
% This is not simple because the layer information has been lost by this point. However,
% infoStr.grayFile should contain the name of the gray file that was used to
% make the flat.mat file. This gray file can be used to identify l1 gNodes.
% If it doesn't exist, the user is prompted to browse for a gray file.
% You can enter this routine with no parameters in which case you'll be asked
% to browse for the input and output files as well
%
% Output files contain an extra field in the infoStr: gloc3d_fixed (==indices of fixed nodes if the fix was successful, 0 otherwise)

if (~exist('inFileName','var'))
    [inFileName,inPathName]=uigetfile('','Locate a flat file');
    inFileName=[inPathName,inFileName];
end

flatFile=load(inFileName);

if (~isfield(flatFile,'gLocs3d'))
    error('No gLocs3d in this file. Can''t go on');
    
end

gLocs3d=flatFile.gLocs3d;

if (~isfield(flatFile,'infoStr'))
    disp('Warning - no infoStr field found.');
    thisInfo='';
else
    thisInfo=flatFile.infoStr;
end

% Check version information on this file to see if it needs to be fixed.
if (isfield(thisInfo,'Version'))
    if (thisInfo.Version>1.5)
        disp('This file is generated by mrFlatMesh > 1.5. It does not need to be fixed.');
        return;    
    end
else
    disp('No version info in file. It may need fixing');
end


% Okay, try to load in the gray file

grayFileName=thisInfo.grayFile;
disp('Attempting to load in the gray file that this flat.mat was made from..');
[nodes, edges, vSize] = readGrayGraph(grayFileName);
if (~nodes & ~edges & ~vSize) % This happens when the file could not be read
    disp(['Failed to read the gray file: ',grayFileName]);
    disp('This may be because you are running on a different platform. Or the gray file may have moved.');
    disp('Try to locate it now....');
    [grayFileName,grayPathName]=uigetfile('','Locate a gray file');
    grayFileName=[grayPathName,grayFileName];
    if (~grayFileName)
        disp('Exiting');
        return;
    else % if a gray file was selected
        % try to load this in..
        [nodes, edges, vSize] = readGrayGraph(grayFileName);
        if (~nodes & ~edges & ~vSize) % This happens when the file could not be read
            error ('Could not read this gray file');
        end % endif bad gray file read
    end % endif gray file selected with uigetfile
end % endif couldn't load the gray file from the infoStr

% We now have nodes, edges, vSize
% nodes(6,:) give the layer...
l1NodeIndices=find(nodes(6,:)==1);
l1NodeCoords=nodes(1:3,l1NodeIndices);
% We now have the coordinates of the l1 nodes in the full gray. We need to know 
% which of these are also in gLocs3d

% Note that this isn't as simple as it sounds: we can subtract [1 1 1] from nodes and then
% do an intersection with gLocs3d but we're also going to catch some l2,3,4 nodes as well
% The only thing  that saves us is that we know that gLocs3d was created by concatenating
% l1,l2,l3,l4 so we know that all the l1 nodes should be in a contiguous block at the
% start. We might get a few errors at the end of the block though.
%
% First add 1 to all of gLocs3d
gLocsInc3d=gLocs3d+1;

% Now use sub2ind to prepare for the intersect
maxDim=max([gLocsInc3d(:);l1NodeCoords(:)]);
gLocsIncInd=sub2ind([maxDim,maxDim,maxDim],gLocsInc3d(:,1),gLocsInc3d(:,2),gLocsInc3d(:,3));
l1NodeInd=sub2ind([maxDim,maxDim,maxDim],l1NodeCoords(1,:),l1NodeCoords(2,:),l1NodeCoords(3,:));

% Find candidate l1nodes
[candidateL1Nodes,indAB,indBA]=intersect(gLocsIncInd,l1NodeInd);

% This returns the nodes that have the same locations in the gray file and the flat gLocs3d
% indAB is an index vector such that candidateL1Nodes=gLocsIncInd(indAB)
% indBA is an index such that candiadateL1Nodes=l1NodeInd(indBA)
% So: indAB should contain a contiguous run from 1:(numL1GNodes) where numL1GNodes is the number of l1GNodes
% it may also contain extra entries where it intersected with other random l2,l3,l4 nodes

% We only want the initial, contiguous indices.
indAB=indAB(:); % Make sure
indABPlus1=shift(indAB,[1,0]); % Circular shift the entries up by 1
diffAB=indAB-indABPlus1; % This will be 1 for the first set of contiguous entries
endofBlock=find(diffAB~=1);
endofBLock=endofBlock(1);

% We're almost done :)
gLocs3d(1:endOfBlock,:)=gLocs3d(1:endOfBlock)+1;
% TOTO Mark field to say that we've fixed it.
% Save out the new flat file

