anchor.diff:  
       if (z == 1.0_dp) then
           do l = 0, lmax
  -            p(l+1) = sqrt( dble(2*l+1) )
  -            dp1(l+1) = sqrt( dble(2*l+1) ) * dble(l) * dble(l+1) / 2.0_dp
  +            p(l+1) = sqrt( dble(2*l+1)) / sqrt(4*pi)
  +            dp1(l+1) = sqrt( dble(2*l+1) ) * dble(l) * dble(l+1) / 2.0_dp &
  +                        / sqrt(4*pi)
           end do
   
       else if (z == -1.0_dp) then
           do l = 0, lmax
  -            p(l+1) = sqrt( dble(2*l+1) ) * dble((-1)**l)
  +            p(l+1) = sqrt( dble(2*l+1)) * dble((-1)**l) / sqrt(4*pi)
               dp1(l+1) = sqrt( dble(2*l+1)) * dble(l) * dble(l+1) &
  -                        * dble((-1)**(l-1)) / 2.0_dp
  -            end do
  +                        * dble((-1)**(l-1)) / 2.0_dp / sqrt(4*pi)
  +        end do
   
       else
           sinsq = (1.0_dp - z**2)
  
positive.code:         end if
  
      else if (size(dp1) < lmax+1) then
          print*, "Error --- PlBar_d1"
          print*, "DP1 must be dimensioned as (LMAX+1) where LMAX is ", lmax
          print*, "Input array is dimensioned ", size(dp1)
          if (present(exitstatus)) then
  
negative.code:         end if
  
      else if (abs(z) > 1.0_dp) then
          print*, "Error --- PlON_d1"
          print*, "ABS(Z) must be less than or equal to 1."
          print*, "Input value is ", z
          if (present(exitstatus)) then
  