function [header, slo, BScans] = octread (path, show)

% This function performs volume OCT data (xxxx.vol) reading. MATLAB syntax of function call is:
% [header, slo, BScans, ILM, RPE] = octread ('path of vol data');
% clc
% close all
% clear all

%% pro načtení toukového souboru - UPRAVIT
adresar = 'C:\Users\Fierka\Documents\BAKALARKA\data\normalni\'; % cesta k datum bez posledniho podadresare 
ID_scan_number='90789' % nazev slozky, kde jsou ulozeny Meta data a zaroven musi byt s stejny korenovy nazev metadat 
addpath ('C:\Users\Fierka\Documents\BAKALARKA\data\normalni\'); % Cesta do podadresáře s programem pro registraci
addpath ([adresar, ID_scan_number '\']); % Cesty k datům OCT+fundus
DatabasePath = [adresar, ID_scan_number '\'];

DatabaseFiles = dir(DatabasePath);

ValidThicknessFile = [];

for i = 1:size(DatabaseFiles,1)
    if regexpi(DatabaseFiles(i).name, 'thickness_');
        ValidThicknessFile = 1;
    end
end


if ValidThicknessFile
    load(['thickness_' ID_scan_number '.mat']); 
else 
    disp ('Neni vytvoren tloustkovy soubor!!!')
end
%% skript
%b UPRAVIT
path=(['C:\Users\Fierka\Documents\BAKALARKA\data\normalni\' ID_scan_number '\' ID_scan_number '.vol']); %uprava pro skript
show=1; %uprava pro skript



% if nargin==1 % slouzi pro funkci 
%     show = 0;
% end

fid = fopen(path);

%% File header read
header.Version = fread( fid, 12, 'int8' );
header.SizeX = fread( fid, 1, 'int32' );
header.NumBScans = fread( fid, 1, 'int32' );
header.SizeZ = fread( fid, 1, 'int32' );
header.ScaleX = fread( fid, 1, 'double' );
header.Distance = fread( fid, 1, 'double' );
header.ScaleZ = fread( fid, 1, 'double' );
header.SizeXSlo = fread( fid, 1, 'int32' );
header.SizeYSlo = fread( fid, 1, 'int32' );
header.ScaleXSlo = fread( fid, 1, 'double' );
header.ScaleYSlo = fread( fid, 1, 'double' );
header.FieldSizeSlo = fread( fid, 1, 'int32' );
header.ScanFocus = fread( fid, 1, 'double' );
header.ScanPosition = fread( fid, 4, 'schar' );
header.ExamTime = fread( fid, 2, 'int32' );
header.ScanPattern = fread( fid, 1, 'int32' );
header.BScanHdrSize = fread( fid, 1, 'int32' );
header.ID = fread( fid, 16, 'int8' );
header.ReferenceID = fread( fid, 16, 'int8' );
header.PID = fread( fid, 1, 'int32' );
header.PatientID = fread( fid, 21, 'int8' );
header.Padding = fread( fid, 3, 'int8' );
header.DOB = fread( fid, 1, 'double' );
header.VID = fread( fid, 1, 'int32' );
header.VisitID = fread( fid, 24, 'int8' );
header.VisitDate = fread( fid, 1, 'double' );
header.Spare = fread( fid, 1840, 'int8' );

disp(['           Version: ' char(header.Version')]); 
disp(['             SizeX: ' num2str(header.SizeX)]);
disp(['         NumBScans: ' num2str(header.NumBScans)]);
disp(['             SizeZ: ' num2str(header.SizeZ)]);
disp(['            ScaleX: ' num2str(header.ScaleX) ' mm']);
disp(['          Distance: ' num2str(header.Distance) ' mm']);
disp(['            ScaleZ: ' num2str(header.ScaleZ) ' mm']);
disp(['          SizeXSlo: ' num2str(header.SizeXSlo)]);
disp(['          SizeYSlo: ' num2str(header.SizeYSlo)]);
disp(['         ScaleXSlo: ' num2str(header.ScaleXSlo) ' mm']);
disp(['         ScaleYSlo: ' num2str(header.ScaleYSlo) ' mm']);
disp(['FieldSizeSlo (FOV): ' num2str(header.FieldSizeSlo) '°']);
disp(['         ScanFocus: ' num2str(header.ScanFocus) ' dpt']);
disp(['      ScanPosition: ' char(header.ScanPosition')]); %OD - right eye; OS - left eye
% disp(['          ExamTime: ' datestr(header.ExamTime(1)/(1e7*60*60*24)+1601*365.25)]);
disp(['       ScanPattern: ' num2str(header.ScanPattern)]);
disp(['      BScanHdrSize: ' num2str(header.BScanHdrSize) ' bytes']);
disp(['                ID: ' char(header.ID)']);
disp(['       ReferenceID: ' char(header.ReferenceID')]);
disp(['               PID: ' num2str(header.PID)]);
disp(['         PatientID: ' char(header.PatientID')]);
disp(['               DOB: ' datestr(header.DOB+693975)]);
disp(['               VID: ' num2str(header.VID)]);
disp(['           VisitID: ' char(header.VisitID')]);
disp(['         VisitDate: ' datestr(header.VisitDate+693975)]);


status = fseek( fid, 2048, -1 );


%% SLO image read
slo = fread( fid, header.SizeXSlo*header.SizeYSlo, 'uint8' );
slo = reshape( slo, header.SizeXSlo, header.SizeYSlo );
slo = imrotate( slo, 90 );
slo=slo/max(max(slo)); % double

figure(1)
imshow(slo,[])
imwrite(slo,'slo.png')
slo = flipud( slo );
figure(1); subplot(1,2,1);
imshow(slo,[]); 
hold on
scrsz = get(0,'ScreenSize');
figure('Position',[1 0 scrsz(3) scrsz(4)-70]);
clf;
if show==1
    figure(2);subplot(1,2,1);
    imshow(slo,[]);
end

%% B-scans read (including image data and scanning position)
status = fseek( fid, 2048+(header.SizeXSlo*header.SizeYSlo), -1 );

BScans=zeros(header.SizeZ,header.SizeX,header.NumBScans);
ILM=zeros(header.NumBScans-1,header.SizeX);
RPE=zeros(header.NumBScans-1,header.SizeX);

if header.NumBScans == 1
    NScans=1;
else
    NScans=header.NumBScans-1;  
end

    for ii = 1:NScans
%         ii
        %Evaluation of BScan position in SLO image
        status = fseek( fid, 16+2048+(header.SizeXSlo*header.SizeYSlo)+(ii*(header.BScanHdrSize+header.SizeX*header.SizeZ*4)), -1 );
        StartX = fread( fid, 1, 'double' );
        StartX_px=round(StartX/(header.SizeXSlo*header.ScaleXSlo)*header.SizeXSlo); %StartX in pixels
        StartY = fread( fid, 1, 'double' );
        StartY_px=round(StartY/(header.SizeYSlo*header.ScaleYSlo)*header.SizeYSlo); %StartY in pixels
        EndX = fread( fid, 1, 'double' );
        EndX_px=round(EndX/(header.SizeXSlo*header.ScaleXSlo)*header.SizeXSlo); %EndX in pixels
        EndY = fread( fid, 1, 'double' );
        EndY_px=round(EndY/(header.SizeYSlo*header.ScaleYSlo)*header.SizeYSlo); %EndY in pixels
        NumSeg = fread( fid, 1, 'int32' );
        subplot(1,2,1);
        if show==1
            figure(2);
            subplot(1,2,1); line([StartX_px EndX_px],[StartY_px EndY_px]);
        end
        
        %BScan reading
        status = fseek( fid, header.BScanHdrSize+2048+(header.SizeXSlo*header.SizeYSlo)+(ii*(header.BScanHdrSize+header.SizeX*header.SizeZ*4)), -1 );
        oct = fread( fid, header.SizeX*header.SizeZ, 'float32' );
        oct = power( reshape( oct, header.SizeX, header.SizeZ ), 0.25 );
        oct = flipud( oct );
        oct = imrotate( oct, -90 );
        BScans(:,:,ii)=oct;
        subplot(1,2,2);
        if show ==1
            figure(2);
            subplot(1,2,2); imshow( BScans(:,:,ii) );
            drawnow
        end
    % Segmented data read
        status = fseek( fid, 256+2048+(header.SizeXSlo*header.SizeYSlo)+(ii*(header.BScanHdrSize+header.SizeX*header.SizeZ*4)), -1 );
        seg = (fread( fid, NumSeg*header.SizeX, 'float' ))';
        ILM(ii,:) = seg(1:header.SizeX);
        RPE(ii,:) = seg(header.SizeX+1:header.SizeX*2);
        pause


    end

fclose( fid );