function [ data , addInfo] = readBrukerRaw(Acqp, varargin)
%function [ data , addInfo ] = readBrukerRaw(Acqp, [path_to_DataFile], ['specified_NRs', your_array], ['specified_Jobs', your_array], ['precision', precision_string])
%IN:
%   Acqp: the struct generated by an object of Experiment or the
%       function readAcqpFile(...)
%
%
%   Optional [] input arguments:
%   path_to_dataFile: string in '' with relative path to fid-data incl. Filename !
%   
%   Additional [] input arguments:
%   specified_XXX: there are 2 cases: 
%       1.) if you are using an standard-fiFile: you can specify a list of NR's to get read, the first NR starts with 1 
%           eg.: <<read_dataFile_ToMatrixStruct(Acqp,
%           path_to_dataFile,'specified_NRs',[2 5 7]);>> -> only NR's 2, 5 and 7 are read
%       2.) if you are using "jobs" in your acquisition process, you can
%           give a list of jobs you want to read, the first job starts with 0. Its the same procedure
%           like the NR's above. You want to read only the fid-file and no
%           job, choose specified_Jobs=[ -1 ];
%   precision: you can define the precision of the read data: you can
%       choose between single precision 'single' or double precision
%       'double', give 'precision' followed by your choose.
%       eg.: read_Raw_ToMatrixStruct(Acqp, path_to_dataFile, 'precision', 'double');
%           or read_Raw_ToMatrixStruct(Acqp, path_to_dataFile, 'precision', 'single');
%       recomanded for big files: single
%       default: double
%
%Out:
%   data: Cellarry, every cell contains a simple 3D Matrix with the Dimensions: (Channel, Scanvalues, NumberOfScans)  
%         if the directory folder contains no job-files: data is a simple 3D Matrix within a Cell with the Dimensions: (Scanvalues, NumberOfScans, Channel)  
%         if the directory folder contains only job-files: {Job0, Job1, Job2, ...}
%         if the directory folder contains fid-file AND job-files: {fidFile, Job1, Job2, ...}



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Copyright (c) 2012
% Bruker BioSpin MRI GmbH
% D-76275 Ettlingen, Germany
%
% All Rights Reserved
%
% $Id: readBrukerRaw.m,v 1.1 2012/09/11 14:22:10 pfre Exp $
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%



%----------------------------------------------------------------
    %% Define default-value if necesseary:
    
    % Set Parameters of InputParser (Matlab-function)
    [varargin, specified_NRs]=bruker_addParamValue(varargin, 'specified_NRs', '@(x) isnumeric(x)', []);
    [varargin, specified_Jobs]=bruker_addParamValue(varargin, 'specified_Jobs', '@(x) isnumeric(x)', []);
    [varargin, precision]=bruker_addParamValue(varargin, 'precision', '@(x) (strcmpi(x,''double'') || strcmpi(x,''single''))', 'double');
    
    if length(varargin) == 1
        path_to_dataFile=varargin{1};
    elseif isempty(varargin) 
        path_to_dataFile='/fid   ';
    else
        warning('MATLAB:bruker_warning', 'Check your input arguments of function readBrukerRaw')
    end
    
%----------------------------------------------------------------   
     %% Check for missing variables in structs:
     cellstruct{1}=Acqp;
     all_here = bruker_requires(cellstruct, {{'Acqp','GO_raw_data_format','BYTORDA','NI','NR','ACQ_size','GO_data_save','GO_block_size', 'AQ_mod'}});
    if isfield(Acqp, 'ACQ_jobs_size')
        all_here = bruker_requires(cellstruct, {{'Acqp','ACQ_jobs','ACQ_jobs_size'}}) & all_here;
    end
    clear cellstruct;
    if ~all_here
        error('Some parameters are missing');
    end

%----------------------------------------------------------------
    %% copy necessary parameters from inputstructs\
    % under this section, the Acqp-struct is not used !
    % makes it easier to find failing parameters for users
    GO_raw_data_format=Acqp.GO_raw_data_format;
    BYTORDA=Acqp.BYTORDA;
    AQ_mod=Acqp.AQ_mod;
    
    % Determining number of selected receive channels
    numSelectedReceivers = bruker_getSelectedReceivers(Acqp);
    
    addInfo.numSelectedReceivers=numSelectedReceivers;
    if ~isempty(specified_NRs)
        addInfo.specified_NRs=specified_NRs;
    end
    if ~isempty(specified_Jobs)
        addInfo.specified_Jobs=specified_Jobs;
    end
    
    %saving the variables to an struct to make function-calls shorter
    paramStruct.NI=Acqp.NI;
    paramStruct.NR=Acqp.NR;
    paramStruct.ACQ_size=Acqp.ACQ_size;
    paramStruct.GO_data_save=Acqp.GO_data_save;
    paramStruct.GO_block_size=Acqp.GO_block_size;
    
    % If Jobs copy variables:
    if isfield(Acqp, 'ACQ_jobs_size') && Acqp.ACQ_jobs_size>0
        jobsExist=true;
        ACQ_jobs=Acqp.ACQ_jobs;
        ACQ_jobs_size=Acqp.ACQ_jobs_size;
        paramStruct.ACQ_jobs=Acqp.ACQ_jobs;
        paramStruct.ACQ_jobs_size=ACQ_jobs_size;
    else
        jobsExist=false;
    end
     
    %----------------------------------------------------------------------
    %% Transform Variables from Acqp-Struct to Matlab readable:

    %transform the Number-Format to Matlab-format (=format) and save the
    %number of bits per value
    switch(GO_raw_data_format)
        case ('GO_32BIT_SGN_INT')
            format='int32';
            bits=32;
        case ('GO_16BIT_SGN_INT')
            format='int16';
            bits=16;
        case ('GO_32BIT_FLOAT')
            format='float32';
            bits=32;
        otherwise
            format='int32';
            disp('Data-Format not correct specified! Set to int32')
            bits=32;
    end

    %transform the machinecode-format to matlab-format (=endian)
    switch(BYTORDA)
        case ('little')
            endian='l';
        case ('big')
            endian='b';
        otherwise 
            endian='l';
            disp('MacineCode-Format not correct specified! Set to little-endian')
    end
    
    % decide if RawFile is complex or real:
    switch AQ_mod
        case ('qf')
            isComplexRaw=false;
        case ('qseq')
            isComplexRaw=true;
        case ('qsim')
            isComplexRaw=true;
        case ('qdig')
            isComplexRaw=true;
        otherwise
            error('The value of parameter AQ_mod is not supported');
    end
        
    %----------------------------------------------------------------------    
    %% Choose the right function(s) for reading:
       
    if ~jobsExist % Parameter doesn's exist
        % 'normal' read of fid-File:
        [ data{1} ] = readfidFile(paramStruct, path_to_dataFile, isComplexRaw, specified_NRs, precision,format,bits, endian,numSelectedReceivers);
        
    elseif (ACQ_jobs(1,1)==0) % if ScanSize entry of First job=0 -> 'normal' read of fid-File:
        [ data{1} ] = readfidFile(paramStruct, path_to_dataFile, isComplexRaw, [], precision,format,bits, endian,numSelectedReceivers);
        
        % special case: first job=normal fidFile, but there are more jobs
        % ACQ_jobs_size>=2 => also jobs in this acquisition
        % if length(specified_Jobs)==1 -> if specified_Jobs(1) is not -1 (-1=read no jobs) => start reading
        % if length(specified_Jobs)==1 => start reading
        if (ACQ_jobs_size>=2) && ( ( length(specified_Jobs)==1 && ~(specified_Jobs(1)==-1 ) || ~length(specified_Jobs)==1 ))
            if isempty(specified_Jobs) % if no jobs specified: choose all jobs (without job0)
                specified_Jobs=[1:ACQ_jobs_size-1];
            else
                specified_Jobs( (specified_Jobs==0) ) = []; % delete the zero-job if an fidFile exist and the user specified an  
            end        
            [ JobsCell ] = readJobFiles(paramStruct, path_to_dataFile, isComplexRaw, specified_Jobs, precision,format,bits, endian,numSelectedReceivers); %=Jobstruct
            data{2:length(JobsCell)} = JobsCell{2:end}; % in this case: JobsCell{1}==[]
            clear JobsCell;
        end
        
        % ~(ACQ_jobs(1,1)==0 => only jobs in this acquisition
        % if length(specified_Jobs)==1 -> if specified_Jobs(1) is not -1 (-1=read no jobs) => start reading
        % if length(specified_Jobs)==1 => start reading
    elseif (~(ACQ_jobs(1,1)==0)) && ( ( length(specified_Jobs)==1 && ~(specified_Jobs(1)==-1 ) || ~length(specified_Jobs)==1 ))
            [ data ] = readJobFiles(paramStruct, path_to_dataFile, isComplexRaw, specified_Jobs, precision,format,bits, endian,numSelectedReceivers);
    else
        error('Your AcqpFile has an unallowed job-description, or you specification was not correct');
    end
    
end

function [ fidFile ] = readfidFile(paramStruct, path_to_dataFile, isComplexRaw, specified_NR, precision,format,bits, endian,numSelectedReceivers)
    % check if specified_NR came with input -> set a boolean
    if(isempty(specified_NR))
        isNR_specified=false;
    else
        isNR_specified=true;
    end
    
    %path_to_dataFile=[path_to_dataFile, '/fid'];
    
    % check the value of precision and change it to a boolean
    if(strcmpi('double',precision))
        single_bool=false;
    elseif(strcmpi('single',precision))
        single_bool=true;
    else
        single_bool=false;
        disp('Your precision-input is not correct! Set to double');
    end      
    clear('input', 'precision');
    
% %----------------------------------------------------------------

    %% Failure Check: is data stored?
    
    % check if acqp-parameter is st correctly
    if(strcmpi(paramStruct.GO_data_save,'no'))
            error('myApp:argChk','You didn''t stored your acqusition! Please change the parameters in ParaVision to save it.');
    end
%----------------------------------------------------------------
    
    %% Transform Variables from paramStruct-Struct to Matlab readable:
    
    % shorter variable-names: 
    NI=paramStruct.NI;
    NR=paramStruct.NR;
    ACQ_size=paramStruct.ACQ_size;
    
    % Calculating number of elements in higher dimensions
    numDataHighDim=prod(paramStruct.ACQ_size(2:end));
    
    % Calculating block size (non-complex)
    if strcmp(paramStruct.GO_block_size,'Standard_KBlock_Format')
        blockSize = ceil(ACQ_size(1)*numSelectedReceivers*(bits/8)/1024)*1024/(bits/8);
    else
        blockSize = ACQ_size(1)*numSelectedReceivers;
    end
    
    
%----------------------------------------------------------------

    %% Read Process:
    % open file       
    try
        fileID = fopen(path_to_dataFile,'r');
    catch
        fileID = -1;
    end
    if fileID == -1
        error('Cannot open parameter file. Problem opening file %s.',path_to_dataFile);
    end
       
    % check if specific NR's are choosen:
    if isNR_specified     
        %% Special case: read only an specified array of NRs:        
        
        %check if maximum user-selected NR > real NR or if it is an bad
        %value
        if (max(specified_NR)>NR || min(specified_NR) < 1 )
            error('Your selected NR is to high or to low. Don''t forget the smallest value is 1. -> function abort !');
        end
        
        % reconfigure NR (for errorcheck and sorting later):
        NR=length(specified_NR); % (don't save to acqp)
        
        %calculate Size of ONE NR, defined as: blockSize * ACQ_sizes-starting-by-dim2 * NI
        size_NR=blockSize*numDataHighDim*NI;
        
        %calculate matrix-dimensions of ONE NR non-komplex:
        num_columns=numDataHighDim*NI; %Spalten
        num_rows=blockSize; %Zeilen
        
        
        
        % Read process:
        %--------------
        
        %Init:
        if single_bool
            fidFile=zeros( num_rows, length(specified_NR)*num_columns, 'single'); % fidFile=temporary variable for fidFile
        else
            fidFile=zeros( num_rows, length(specified_NR)*num_columns); % fidFile=temporary variable for fidFile
        end
        
        % read th specified NR-Data:        
        for i=1:length(specified_NR)
            % set position in file:
            fseek(fileID, ( specified_NR(i)-1 )*size_NR*bits/8 , 'bof');
            
            %read only one NR with single or double precision:
            if(single_bool)
                Xzw=single(fread(fileID, [num_rows, num_columns], [format, '=>single'], 0 , endian));
            else
                Xzw=fread(fileID, [num_rows, num_columns], format, 0 , endian);
            end         
            % write to the output-variable:
            fidFile(:,(i-1)*num_columns+1:i*num_columns)=Xzw(1:end,:);
        end        
    else 
        %% normal-case: read complete fid-file
        
        %read File to Variable fidFile with fread() and make it single-precission with single():
        %  Attention: Matlab ACQ_size(1) = ACQ_size(0)
        if(single_bool)
            fidFile=( fread(fileID, [blockSize, numDataHighDim*NI*NR], [format, '=>single'], 0, endian) );
        else
            fidFile=fread(fileID, [blockSize, numDataHighDim*NI*NR], format, 0, endian);
        end    
    end
    
    % file close
    fclose(fileID);
    
    %% short errorcheck and cut the zeros
    fidFileSize=numel(fidFile);
    if ~(fidFileSize==blockSize*numDataHighDim*NI*NR)
        error('Size of fid file does not match parameters.');
    end
        
    fidFile=reshape(fidFile,blockSize,numDataHighDim*NI*NR); % in most cases unnecessary
   
    % for faster execution: minimize usage of the permute-command
    if blockSize ~= ACQ_size(1)*numSelectedReceivers %remove zero-lines
        fidFile=permute(fidFile,[2,1]); %for faster memory access during 2 following operations
        fidFile=fidFile(:,1:ACQ_size(1)*numSelectedReceivers);
        % select channels for new dimension:  
        fidFile=reshape(fidFile,[numDataHighDim*NI*NR, ACQ_size(1) ,numSelectedReceivers]); 
           
        % resort dimensions
        fidFile=permute(fidFile,[3,2,1]);
    else 
        fidFile=reshape(fidFile,[ACQ_size(1), numSelectedReceivers, numDataHighDim*NI*NR]);
        fidFile=permute(fidFile,[2 1 3]);
    end
    
    
    
    % Save to output variable:
    if isComplexRaw
        % convert to complex:
        fidFile=complex(fidFile(:,1:2:end,:,:), fidFile(:,2:2:end,:,:));
    % else: don't convert, only save
    end
    
end

function [ JobsCell ] = readJobFiles(paramStruct, path_to_dataFile, isComplexRaw, specified_job, precision,format,bits, endian,numSelectedReceivers)
% specified_jobs or empty entry:
    if(isempty(specified_job))
        jobs=[0:paramStruct.ACQ_jobs_size-1];
    else
        jobs=specified_job;
    end
    clear specified_job;
    
    if max(jobs)>paramStruct.ACQ_jobs_size-1 || min(jobs)<0
        error('Your specified jabs are not correct !')
    end
    
    % check the value of precision and change it to a boolean
    if(strcmpi('double',precision))
        single_bool=false;
    elseif(strcmpi('single',precision))
        single_bool=true;
    else
        single_bool=false;
        disp('Your precision-input is not correct! Set to double');
    end
    
    clear('input', 'precision');
    
% %----------------------------------------------------------------

    %% Failure Check: is data stored?
    
    % check if acqp-parameter is st correctly
    if(strcmpi(paramStruct.GO_data_save,'no'))
            error('myApp:argChk','You didn''t stored your acqusition! Please change the parameters in ParaVision to save it.');
    end    
 
%----------------------------------------------------------------
    %% read jobs:
    for i=jobs
        % Modify path: before '../yourpath/fid'
        pos=strfind(path_to_dataFile, '/'); % last '/' is end of directory description
        path_to_dataFile=path_to_dataFile(1:pos(end)); % -> after: '../yourpath/'
        temppath_to_dataFile=[path_to_dataFile, 'rawdata.job', num2str(i)];
        
        %% Read Process:
        % open file       
        try
            fileID = fopen(temppath_to_dataFile,'r');
        catch
            fileID = -1;
        end
        if fileID == -1
            error('Cannot open parameter file. Problem opening file %s.',temppath_to_dataFile);
        end


        %read File to Variable X with fread() and make it single-precission with single():
        %  Attention: Matlab ACQ_size(1) = ACQ_size(0)
        if(single_bool)
            X=( fread(fileID, [numSelectedReceivers*paramStruct.ACQ_jobs(1,i+1), inf], [format, '=>single'], 0, endian) );
        else
            X=fread(fileID, [numSelectedReceivers*paramStruct.ACQ_jobs(1,i+1), inf], format, 0, endian);
        end
        
        dim1=numel(X) / (paramStruct.ACQ_jobs(1,i+1)*numSelectedReceivers);
        X=reshape(X,[numSelectedReceivers, paramStruct.ACQ_jobs(1,i+1),dim1]);
        
%         % select channels for new dimension:
%         dim1=numel(X) / (paramStruct.ACQ_jobs(1,i+1)*numSelectedReceivers); % = equal to numDataHighDim*NI*NR
%         X=reshape(X,[dim1, paramStruct.ACQ_jobs(1,i+1), numSelectedReceivers]); 
% 
%         % resort dimensions
%         X=permute(X,[2,1,3]);
        
        
        % Save to output variable:
        if isComplexRaw
            % convert to complex:
            X=complex(X(:,1:2:end,:), X(:,2:2:end,:));
        % else: don't convert, only save
        end

        % file close
        fclose(fileID);

        % Add to struct:
        
        JobsCell{i+1}=X;
    end

end
