function dataTable = importJeti(filename, selectedDatasets)
%IMPORTJETI Import JETI SPECBOS irradiance data (CSV or TXT format)
%  DATATABLE = IMPORTJETI(FILENAME) reads all data from the text file.
%
%  DATATABLE = IMPORTJETI(FILENAME, SELECTEDDATASETS) reads only the 
%  Wavelength and the specific dataset indices you request.
%  Example: table = importJeti("20260421_1500_Jeti.csv", [1, 5, 10]);
%
%  Supports:
%  1. CSV: Auto-detects comma (,) or semicolon (;) delimiters, handles
%     various Date/Time formats, and flexible header lengths.
%  2. isd: txt-style format with [sections] and tab-separated data

if nargin < 2
    selectedDatasets = []; % Empty means import all
end

%% Determine file type and route accordingly
[~, ~, ext] = fileparts(filename);

if strcmpi(ext, '.csv')
    dataTable = importCSV(filename, selectedDatasets);
elseif strcmpi(ext, '.isd')
    dataTable = importTXT(filename);
else
    error('Unsupported file format: %s. Use .csv or .isd', ext);
end

end

%% ========== CSV IMPORT FUNCTION ==========
function dataTable = importCSV(filename, selectedDatasets)

fid = fopen(filename, 'r');
if fid == -1
    error('Cannot open file: %s', filename);
end

% 1. Auto-detect Delimiter
fileDelimiter = ','; % default
while ~feof(fid)
    thisLine = fgetl(fid);
    if contains(thisLine, ';')
        fileDelimiter = ';';
        break;
    elseif contains(thisLine, ',')
        fileDelimiter = ',';
        break;
    end
end
frewind(fid); % Reset back to beginning of file

metadata = struct();
metadata.numDatasets    = NaN;
metadata.numPoints      = NaN;
metadata.Date           = {};
metadata.Time           = {};

lineIdx            = 0;
dataHeaderLineIdx  = NaN;
nHeaderCols        = NaN;

% 2. Parse Metadata and find Wavelength header
while true
    thisLine = fgetl(fid);
    if ~ischar(thisLine)
        break
    end
    lineIdx = lineIdx + 1;

    parts = strsplit(thisLine, fileDelimiter);
    if isempty(parts)
        continue
    end

    % Detect start of data header line (usually 'Wavelength [nm]')
    if contains(parts{1}, 'Wavelength', 'IgnoreCase', true) && ...
       ~contains(parts{1}, 'Start Wavelength') && ...
       ~contains(parts{1}, 'End Wavelength')
        dataHeaderLineIdx = lineIdx;
        nHeaderCols = numel(parts);
        break
    end

    key = strtrim(parts{1});

    switch key
        case 'Number of datasets'
            if numel(parts) >= 2
                metadata.numDatasets = str2double(strtrim(parts{2}));
            end
        case 'Number of points'
            if numel(parts) >= 2
                metadata.numPoints = str2double(strtrim(parts{2}));
            end
        case 'Date'
            if numel(parts) >= 2
                % Remove empty trailing cells
                metadata.Date = parts(2:end);
                metadata.Date = metadata.Date(~cellfun('isempty', metadata.Date));
            end
        case 'Time'
            if numel(parts) >= 2
                metadata.Time = parts(2:end);
                metadata.Time = metadata.Time(~cellfun('isempty', metadata.Time));
            end
    end
end

if isnan(dataHeaderLineIdx)
    fclose(fid);
    error('Could not find Wavelength data header line in CSV.');
end

firstDataLineInFile = dataHeaderLineIdx + 1;
fclose(fid);

% 3. Determine column structures
if isnan(metadata.numDatasets)
    % Infer from header columns (Subtract 1 for Wavelength column)
    metadata.numDatasets = nHeaderCols - 1;
end

numCols = 1 + metadata.numDatasets;

% Build dataLines bounds
if isnan(metadata.numPoints)
    dataLines = [firstDataLineInFile, Inf];
else
    lastDataLineInFile = firstDataLineInFile + metadata.numPoints - 1;
    dataLines = [firstDataLineInFile, lastDataLineInFile];
end

% Set up Import Options
opts = delimitedTextImportOptions("NumVariables", numCols);
opts.DataLines = dataLines;
opts.Delimiter = fileDelimiter;

% 4. Build Variable Names
varNames = cell(1, numCols);
varNames{1} = "Wavelength [nm]";

if ~isempty(metadata.Date) && ~isempty(metadata.Time)
    nMeta = min([numel(metadata.Date), numel(metadata.Time), metadata.numDatasets]);
    for k = 1:nMeta
        dstr = strtrim(metadata.Date{k});
        tstr = strtrim(metadata.Time{k});
        
        timestr = parseJetiDateTime(dstr, tstr);
        if isempty(timestr)
            timestr = sprintf('dataset%02d', k);
        end

        varNames{1 + k} = sprintf("SpectralIrradiance_%s [W/(sqm*nm)]", timestr);
    end

    for k = nMeta+1:metadata.numDatasets
        varNames{1 + k} = sprintf("SpectralIrradiance_%02d [W/(sqm*nm)]", k);
    end
else
    for k = 1:metadata.numDatasets
        varNames{1 + k} = sprintf("SpectralIrradiance_%02d [W/(sqm*nm)]", k);
    end
end

opts.VariableNames = string(varNames);
opts.VariableTypes = repmat("double", 1, numCols);
opts.ExtraColumnsRule = "ignore";
opts.EmptyLineRule   = "read";
opts.VariableNamingRule = "preserve";

% Handle Specific Selected Columns
if ~isempty(selectedDatasets)
    colsToKeep = [1, selectedDatasets + 1]; 
    % Prevent out-of-bounds error
    colsToKeep = colsToKeep(colsToKeep <= numCols);
    opts.SelectedVariableNames = opts.VariableNames(colsToKeep);
end

% 5. Read the Table
dataTable = readtable(filename, opts);

end

%% ========== DATETIME PARSER HELPER ==========
function timestr = parseJetiDateTime(dstr, tstr)
% Tries several Date/Time formats to safely parse standard 24h, 12h(am/pm), 
% and regional variations (dd.MM.yyyy vs MM/dd/yyyy).
    combinedStr = [dstr ' ' tstr];
    
    formatsToTry = {
        'MM/dd/yyyy hh:mm:ssa', ... % US 12-hour (03:05:17pm)
        'dd.MM.yyyy HH:mm:ss',  ... % EU 24-hour (12.12.2025 15:54:00)
        'MM/dd/yyyy HH:mm:ss',  ... % US 24-hour
        'dd/MM/yyyy HH:mm:ss',  ... % EU slash 24-hour
        'yyyy-MM-dd HH:mm:ss'       % ISO 8601
    };

    timestr = '';
    for i = 1:numel(formatsToTry)
        try
            dt = datetime(combinedStr, 'InputFormat', formatsToTry{i});
            if ~isnat(dt)
                timestr = datestr(dt, 'yyyymmdd_HHMMss');
                return;
            end
        catch
            % Format mismatched, try the next one
            continue;
        end
    end
end

%% ========== TXT IMPORT FUNCTION ==========
function dataTable = importTXT(filename)

fid = fopen(filename, 'r');
if fid == -1
    error('Cannot open file: %s', filename);
end

metadata = struct();
metadata.Name           = '';
metadata.NumberOfDataX  = NaN;
metadata.NumberOfDataY  = NaN;
metadata.XMin           = NaN;
metadata.XMax           = NaN;
metadata.XUnit          = '';
metadata.YUnit          = '';

dataStartLineIdx = NaN;
lineIdx = 0;

% Parse header sections
while true
    thisLine = fgetl(fid);
    if ~ischar(thisLine)
        break
    end
    lineIdx = lineIdx + 1;
    
    if strcmp(strtrim(thisLine), 'Data')
        dataStartLineIdx = lineIdx + 1;
        break
    end
    
    if contains(thisLine, '=')
        parts = strsplit(thisLine, '=');
        if numel(parts) == 2
            key   = strtrim(parts{1});
            value = strtrim(parts{2});
            
            switch key
                case 'Name'
                    metadata.Name = value;
                case 'NumberOfDataX'
                    metadata.NumberOfDataX = str2double(value);
                case 'NumberOfDataY'
                    metadata.NumberOfDataY = str2double(value);
                case 'X-MIN'
                    metadata.XMin = str2double(value);
                case 'X-MAX'
                    metadata.XMax = str2double(value);
                case 'X-Unit'
                    metadata.XUnit = value;
                case 'Y-Unit'
                    metadata.YUnit = value;
            end
        end
    end
end

if isnan(dataStartLineIdx)
    fclose(fid);
    error('Could not find data section in TXT file.');
end

fclose(fid);

opts = delimitedTextImportOptions("NumVariables", 2);
opts.DataLines = [dataStartLineIdx, Inf];
opts.Delimiter = '\t';
opts.VariableNamingRule = "preserve";

if ~isempty(metadata.XUnit) && ~isempty(metadata.YUnit)
    varNames = string({
        sprintf("Wavelength [%s]", metadata.XUnit)
        sprintf("SpectralIrradiance_%s [%s]", metadata.Name, metadata.YUnit)
    });
else
    varNames = string({
        sprintf("Wavelength_%s", metadata.Name)
        sprintf("SpectralIrradiance_%s", metadata.Name)
    });
end

opts.VariableNames = varNames;
opts.VariableTypes = ["double", "double"];

dataTable = readtable(filename, opts);
end