%--------------------------------------------------------------------------
% Typical Mapped Spectral Accelerations (ASCE/SEI 7-22 §11.4.1)
%
%   Seismicity Zone     Example Locations                       Typical Ss (g)      Typical S1 (g)
%   -----------------   --------------------------------------   ----------------   ----------------
%   Very low            Florida, northern Midwest                0.05 – 0.15         0.02 – 0.05
%   Low–moderate        Huntsville AL, Dallas TX, East Coast     0.25 – 0.60         0.08 – 0.20
%   Moderate–high       St Louis MO (New Madrid), Charleston SC  0.60 – 1.0          0.20 – 0.40
%   High                Salt Lake City UT, Seattle WA            1.0 – 1.5           0.4 – 0.6
%   Very high           Los Angeles CA, San Francisco CA         1.5 – 2.5 (~max 2.75)  0.6 – 1.0 (~max 1.25)
%   Extreme             Anchorage AK, southern Alaska            up to 3.0           up to 1.5
%
%   Notes:
%     • Ss = mapped spectral acceleration at short periods (0.2 s)
%     • S1 = mapped spectral acceleration at 1-second period
%     • Values are approximate, for general reference only.
%       Obtain exact values from the ASCE 7 Hazard Tool or USGS Seismic Design Maps.
%
%--------------------------------------------------------------------------
% Approximate Long-Period Transition Period, TL (s)
%   per ASCE/SEI 7-22 Table 11.4-3 and typical regional trends
%
%   Site Class    Typical Soil Profile Description           TL (s)
%   -----------   ----------------------------------------   ------
%       A         Hard rock                                  8.0
%       B         Rock                                       8.0
%       C         Very dense soil and soft rock              6.0
%       D         Stiff soil profile                         4.0
%       E         Soft clay soil                             2.0
%
%   Notes:
%     • Actual TL values vary by geographic region per ASCE 7-22
%       long-period transition maps (Figures 22-17 – 22-20).
%     • These values are suitable as default or fallback estimates
%       when a mapped TL is not available.
%--------------------------------------------------------------------------


function [T,fn,Sa] = Sa_horizontal_from_Ss(Ss, S1, siteClass, TL)
% Sa_horizontal_from_Ss  Horizontal two-period design response spectrum
%   per ASCE/SEI 7-22 §11.4.5.2 using (Ss, S1, siteClass) inputs.
%
%   Sa = Sa_horizontal_from_Ss(Ss, S1, siteClass, TL)
%
%   Inputs:
%     Ss        = mapped spectral acceleration at 0.2 s (g)  [ASCE/USGS maps]
%     S1        = mapped spectral acceleration at 1.0 s (g)  [ASCE/USGS maps]
%     siteClass = 'A'..'E' (Site Class F requires site-specific study → error)
%     TL        = long-period transition (s) [ASCE long-period map]
%
%   Output:
%     Sa  = spectral acceleration (g) evaluated on a 1/12-octave T grid
%
%   Notes:
%     - Fa(Ss,siteClass) from Table 11.4-1; Fv(S1,siteClass) from Table 11.4-2.
%     - Linear interpolation within tables; clamped to end columns outside ranges.
%     - SDS = (2/3)*Fa*Ss; SD1 = (2/3)*Fv*S1; 5% critical damping convention.
%     - Uses external plotting helpers: plot_loglog_function_h2, plot_loglin_function_h2.
%

    % --- Validate inputs ---
    validateattributes(Ss, {'numeric'}, {'scalar','real','finite','positive'});
    validateattributes(S1, {'numeric'}, {'scalar','real','finite','positive'});
    validateattributes(TL, {'numeric'}, {'scalar','real','finite','positive'});
    if ~(ischar(siteClass) || isstring(siteClass))
        error('siteClass must be a letter A..E'); 
    end
    sc = upper(string(strtrim(siteClass)));
    if sc == "F"
        error('Site Class F requires site-specific procedures per ASCE 7-22.');
    end

    % --- Lookup Fa, Fv via ASCE 7-22 Tables (interpolated) ---
    Fa = fa_table_interp(Ss, sc);
    Fv = fv_table_interp(S1, sc);

    % --- Compute SDS, SD1 ---
    SMS = Fa * Ss;
    SM1 = Fv * S1;
    SDS = (2/3) * SMS;
    SD1 = (2/3) * SM1;

    % --- Title string ---
    t_string = sprintf('Horizontal Design Response Spectrum 5%% Damping \n Site Class %s (§11.4.5.2)', sc);

    % --- 1/12-octave frequency grid ---
    f_min = 0.1; f_max = 100; ratio = 2^(1/12);
    frequencies = f_min;
    f = f_min;
    while f * ratio <= f_max
        f = f * ratio;
        frequencies(end+1,1) = f; %#ok<SAGROW>
    end
    if frequencies(end) < f_max
        frequencies(end+1,1) = f_max;
    end
    T = sortrows(1 ./ frequencies);

    % --- Characteristic periods (§11.4.5.2) ---
    T0 = 0.2 * (SD1 / SDS);    % Eq. 11.4-5
    TS =       (SD1 / SDS);    % (definition)

    % --- Piecewise spectrum ---
    Sa = zeros(size(T));
    idx1 = (T <  T0);              Sa(idx1) = SDS * (0.4 + 0.6 * (T(idx1) ./ T0));
    idx2 = (T >= T0) & (T <= TS);  Sa(idx2) = SDS;
    idx3 = (T >  TS) & (T <= TL);  Sa(idx3) = SD1 ./ T(idx3);
    idx4 = (T >  TL);              Sa(idx4) = (SD1 ./ TL) .* (TL ./ T(idx4)).^2;

    % === Plotting (uses your helpers) ===
    fig_num = 1;
    y_label = 'Spectral acceleration, Sa(T) (g)';

    % Period-domain
    ppp = [T Sa];
    x_label = 'Period T (sec)';
    [fig_num,~] = plot_loglog_function_h2(fig_num, x_label, y_label, t_string, ppp, min(T), max(T));
    [fig_num,~] = plot_loglin_function_h2(fig_num, x_label, y_label, t_string, ppp, min(T), max(T));

    % Frequency-domain
    f = 1 ./ T;
    ppp = flipud([f Sa]);  % increasing freq to the right
    x_label = 'Natural Frequency fn (Hz)';
    [fig_num,~] = plot_loglog_function_h2(fig_num, x_label, y_label, t_string, ppp, min(f), max(f));
    [fig_num,~] = plot_loglin_function_h2(fig_num, x_label, y_label, t_string, ppp, min(f), max(f));

    % --- Optional: print summary ---
    fprintf('\nSite Class %s\n', sc);
    fprintf('  Fa  = %.3f   Fv  = %.3f\n', Fa, Fv);
    fprintf('  SDS = %.3f g SD1 = %.3f g   TL = %.2f s\n\n', SDS, SD1, TL);

    fn=f;
end

%==================== Helpers ====================%

function Fa = fa_table_interp(Ss, sc)
    % ASCE 7-22 Table 11.4-1 (representative values). Linear interpolation.
    Ss_bp = [0.25 0.50 0.75 1.00 1.25 1.50];  % g
    switch sc
        case "A", vals = [0.8 0.8 0.8 0.8 0.8 0.8];
        case "B", vals = [0.8 0.8 0.8 0.8 0.8 0.8];
        case "C", vals = [1.2 1.2 1.1 1.0 1.0 1.0];
        case "D", vals = [1.6 1.4 1.2 1.1 1.0 1.0];
        case "E", vals = [2.5 2.1 1.7 1.3 1.2 1.2];
        otherwise, error('Unsupported site class "%s" for Fa.', sc);
    end
    Fa = interp_clamped(Ss_bp, vals, Ss);
end

function Fv = fv_table_interp(S1, sc)
    % ASCE 7-22 Table 11.4-2 (representative values). Linear interpolation.
    S1_bp = [0.10 0.20 0.30 0.40 0.50 0.60];  % g
    switch sc
        case "A", vals = [0.8 0.8 0.8 0.8 0.8 0.8];
        case "B", vals = [0.8 0.8 0.8 0.8 0.8 0.8];
        case "C", vals = [1.7 1.6 1.5 1.4 1.3 1.3];
        case "D", vals = [2.4 2.2 2.0 1.8 1.7 1.6];
        case "E", vals = [3.5 3.2 2.8 2.4 2.4 2.4];
        otherwise, error('Unsupported site class "%s" for Fv.', sc);
    end
    Fv = interp_clamped(S1_bp, vals, S1);
end

function y = interp_clamped(x_bp, y_bp, xq)
    % Linear interpolation with clamping to end values (ASCE tables convention).
    if xq <= x_bp(1), y = y_bp(1); return; end
    if xq >= x_bp(end), y = y_bp(end); return; end
    y = interp1(x_bp, y_bp, xq, 'linear');
end
