function [cardiac_artifact_shape, peak_offset, var_ratio] = find_cardiac_artifact_shape( vector, r_peak_list, samplerate )
%    vector = nobaseline_ecg;
%    r_peak_list = list_r_peaks;



    PC_COUNT = 1;

    MAX_SEGMENTS_COUNT = 180;

    SECONDS_BEFORE_R = 0.15;
    SECONDS_AFTER_R = 0.5;


    if ~isvector(vector) || ~isnumeric(vector)
        error('cvp:find_cardiac_artifact_shape:vector', 'vector is not a vector');
    end

    if ~isvector(r_peak_list) || ~isnumeric(r_peak_list)
        error('cvp:find_cardiac_artifact_shape:r_peak_list_not_vector', 'r_peak_list is not a numeric vector');
    end

    if (min(r_peak_list) < 1 || max(r_peak_list) > length(vector))
        error('cvp:find_cardiac_artifact_shape:invalid_r_peak_list', 'r_peak_list contains invalid values');
    end

    if (length(r_peak_list) < 2)
        error('cvp:find_cardiac_artifact_shape:not_enough_r_peaks', 'Not enough R-Peaks Provided');
    end

    if ~isscalar(samplerate) || ~isnumeric(samplerate) || (fix(samplerate) ~= samplerate)
        error('cvp:find_cardiac_artifact_shape:samplerate_not_integer', 'samplerate must be an integer');
    end


    row_time = SECONDS_BEFORE_R + SECONDS_AFTER_R;
    hd = create_filter_cheby(8, (1/(row_time*2)) ./ (samplerate / 2), 'high');
    vector = filtfilthd(hd, vector);


    segment_start_offset = floor(SECONDS_BEFORE_R * samplerate);
    peak_offset = segment_start_offset;
    segment_end_offset = floor(SECONDS_AFTER_R * samplerate);

    verified_r_peak_list = r_peak_list(r_peak_list > segment_start_offset & r_peak_list < (length(vector) - segment_end_offset));

    matrix_row_count = min(length(verified_r_peak_list), MAX_SEGMENTS_COUNT);
    r_peak_offset = unique(round(linspace(1,length(verified_r_peak_list), matrix_row_count)));
    matrix_row_count=length(r_peak_offset);


    segment_start = reshape(verified_r_peak_list - segment_start_offset, [], 1); %Make sure it's a Mx1 Vector

    segment_end = reshape(verified_r_peak_list + segment_end_offset, [], 1); %Make sure it's a Mx1 Vector
    segment_length = segment_end(1) - segment_start(1) + 1;

    ecg_mat_ind = segment_start(r_peak_offset) * ones(1,segment_length) + ones(1,matrix_row_count)' * (0:(segment_length-1));
    ecg_mat = vector(ecg_mat_ind);
 

    %TODO: Optimize
    [U,S,V] = svd(ecg_mat, 1);
    %Calc: how much Variance is described by selected PC
    pc_list = diag(S);
    var_ratio = sum(pc_list(1:PC_COUNT)) / sum(pc_list);

    %Avoid calculation if not needed anyways
    if (var_ratio > 0.2)
        S_reconstruct = zeros(size(S));
        S_reconstruct(1,1) = S(1,1);
        cardiac_artifact_shape_mtrx = U * S_reconstruct * V';
        cardiac_artifact_shape = cardiac_artifact_shape_mtrx(1,:);
    %    cardiac_artifact_shape = (U(1,:) .* S(1,1)) * V';
    else
        cardiac_artifact_shape=zeros(1,segment_length);
    end

%% DEBUG TEST TEST TEST
%	var_ratio = 0;
%	cardiac_artifact_shape=zeros(1,segment_length);
end
