function h = hFromT(T,p)
  global p_i;
  global hfTP;
  global dT;
  
    fT = floor(T-273.16+1);
    cT = ceil(T-273.16+1);
    [FP,fP] = min((abs(p - p_i))');
    fP = fP';
    cP = fP + 1;
    
    hf = diag(hfTP(fT,fP)) + (T-273.16+1 - fT)./(cT - fT + 1e-20).*(diag(hfTP(cT,fP)) - diag(hfTP(fT,fP)));
    hc = diag(hfTP(fT,cP)) + (T-273.16+1 - fT)./(cT - fT + 1e-20).*(diag(hfTP(cT,cP)) - diag(hfTP(fT,cP)));
    h = hf + (p - p_i(fP)')./(p_i(cP)' - p_i(fP)' + 1e-20).*(hc - hf);

endfunction
