function [Q_mm_h, meta] = readObservedQ(par)
% READOBSERVEDQ 

    fileInput = getObsParam(par, 'fileInput', 'measuredDischarge.xlsx');
    spdFile = getObsParam(par, 'spdFile', 144);
    units = lower(string(getObsParam(par, 'units', 'm3s')));

    if ~isfile(fileInput)
        error('readObservedQ: file non trovato: %s', fileInput);
    end

    obs = readtable(fileInput);
    qCol = detectQColumn(obs);
    Q_raw = obs.(qCol)(:);

    if units ~= "m3s"
        error('readObservedQ: units "%s" non supportato (usa ''m3s'').', units);
    end

    spdModel = round(24 / par.domainT.dt);
    area_km2 = par.domainL.A;
    ratio = round(spdFile / spdModel);

    if ratio >= 2 && ismember('date', obs.Properties.VariableNames)
        aggTbl = aggregateDischargeTable(obs, qCol, area_km2);
        Q_mm_h = aggTbl.Q_mm_h;
    else
        Q_mm_h = aggregateToModelDt(Q_raw, spdFile, spdModel, 'discharge', area_km2);
    end

    if isfield(par, 'domainT')
        Nt = min(par.domainT.Nt, numel(Q_mm_h));
        Q_mm_h = Q_mm_h(1:Nt);
    end

    meta = struct();
    meta.fileInput = fileInput;
    meta.unit = 'mm/h';
    meta.inputUnit = 'm3/s';
    meta.nSteps = numel(Q_mm_h);
    meta.nDays = meta.nSteps / spdModel;
    meta.spdFile = spdFile;
    meta.spdModel = spdModel;
    meta.isSubdaily = spdFile > spdModel;
    meta.area_km2 = area_km2;

    % fprintf('readObservedQ: %s -> %d passi orari [mm/h] (A = %.2f km^2)\n', ...
    %     fileInput, meta.nSteps, area_km2);
end

function qCol = detectQColumn(obs)
    cols = obs.Properties.VariableNames;
    if ismember('Q', cols)
        qCol = 'Q';
    elseif ismember('Q_obs', cols)
        qCol = 'Q_obs';
    else
        error('readObservedQ: colonna Q o Q_obs non trovata.');
    end
end

function v = getObsParam(par, name, default)
    if isfield(par, 'obs') && isfield(par.obs, name) && ~isempty(par.obs.(name))
        v = par.obs.(name);
    else
        v = default;
    end
end
