subroutine xpxpr_qpar(r0,rr,fc,hm,ym,invqp,f,X,PXPR)
c	s/r to calculate the values of X and PXPR for a single
c	quasi-parabolic layer (invqp = -1 for normal QP, +1 for inverted)
c	written by leo mcnamara 1-jun-88
c-----	last modified 24-oct-90
c
	implicit real*8 (a-h,o-z)
c
	inv = - invqp
	rm = r0 + hm
	rb = rm - ym
c
	x = 0.d0
	pxpr = 0.d0
c
c	the parameter X = (FN/F)**2
	z = (rm - rr) * rb / ym / rr
c	loss of significance will  cause fn2 < 0
	if(dabs(z - 1.0d0).lt.1.d-15) return
c	allow for small errors near the base of the layer
	if(abs(rb-rr).le.0.01d0.and.z.gt.1.d0) z = 1.
	fn2 = fc * fc * (1.d0 - INV * z * z)
	if(fn2.eq.0.d0) return
	if(fn2.lt.0.d0) print*,' XPXPR_QPAR  -- fc,hm,hh,hb,ym,z ',
     a  fc,rm-r0,rr-r0,rb-r0,ym,z
	fn = sqrt (fn2)
	x = fn2 / f / f
	if(x.gt.1.0d0) print*,' rr,fn,f,x',rr,fn,f,x
ccc	if(x.gt.1.0d0) stop
c
c	partial derivative of X wrt r
	pxpfn = 2.d0 * fn / f / f
	pfnpz =  - INV * fc * fc * z / fn
	pzpr = - rb * rm / ym / rr / rr
	pxpr = pxpfn * pfnpz * pzpr
c
c	we always need the positive value
	pxpr = abs(pxpr)
c
	return
	end