function [a_f,b_f,c_f,a_r, b_r, c_r] = coeff_quad_fun(theta_f,theta_r,theta_m,p_m,p_fm,p_rm,p_f,p_r)
    
    % a*theta^2+b*theta+c*theta
    theta_fm= (theta_f+theta_m)/2;
    theta_rm= (theta_r+theta_m)/2;
    % p_f=0; % sigma(theta_f)
    % p_r=0; % sigma(theta_r)

    % parabola decrescente
    M_f= [theta_m^2 theta_m 1;
        theta_fm^2 theta_fm 1;
        theta_f^2 theta_f 1];
    p_vec_f = [p_m; p_fm; p_f];
    x_f= M_f\p_vec_f;

    a_f= x_f(1);
    b_f= x_f(2);
    c_f= x_f(3);

    % parabola crescente
    M_r= [theta_m^2 theta_m 1;
    theta_rm^2 theta_rm 1;
    theta_r^2 theta_r 1];
    p_vec_r = [p_m; p_rm; p_r];
    x_r= M_r\p_vec_r;

    a_r= x_r(1);
    b_r= x_r(2);
    c_r= x_r(3);


end