!
!     quadratic interpolation subroutine , output is alpgt-->interpolated lnP
!
       SUBROUTINE quadInterpLnPr(h1,h2,h3,alp1,alp2,alp3,h,alpgt)
 
       implicit none
       REAL::h1,h2,h3,alp1,alp2,alp3,h,alpgt
       REAL::h3m2,h3m1,a,b,c,d,dlt
       logical debug

         dlt=(alp1-alp3)*(alp2-alp1)*(alp3-alp2)
         h3m2=h3-h2
         h3m1=h3-h1
         b=(h3m2*(alp3*alp3-alp1*alp1)   &
            -h3m1*(alp3*alp3-alp2*alp2))/dlt
         c=(h3m1*(alp3-alp2)-h3m2*(alp3-alp1))/dlt
         d=h1-c*alp1*alp1-b*alp1-h

         if((b*b-4*c*d).ge.0) then
           alpgt=(-b-sqrt(b*b-4*c*d))/2/c
         else
           a=(h1-h2)/(alp1-alp2)
           b=h1-a*alp1
           alpgt=(h-b)/a
         endif

         alpgt = alp1*(h-h2)*(h-h3)/((h1-h2)*(h1-h3)) &
               + alp2*(h-h1)*(h-h3)/((h2-h1)*(h2-h3)) &
               + alp3*(h-h1)*(h-h2)/((h3-h1)*(h3-h2))

      END SUBROUTINE quadInterpLnPr