anchor.diff:  
           ISSUCTION = .FALSE.
           IF (Cf(2) .GT. 0.) THEN
  -            CALL wrap_qk61(omega,f_int2,a,b,answer,abserr,resabs,resasc)   ! wrapper for qk61 routine from slatec library
  +            CALL solve_qk61(omega,f_int2,a,b,answer,abserr,resabs,resasc)   ! wrapper for qk61 routine from slatec library
               Spectrum = D/(4.*pi*R**2.)*answer
               SPL_press = 10.*log10(Spectrum*DBARH/2.e-5/2.e-5)
               SPLP(i_omega) = SPL_press + 10.*log10(band_width)
  
positive.code:     ErrStat = ErrID_None
      ErrMsgn = ""
      ! Body of TNO 
      epsabs = 1e-10     !absolute accuracy
      epsrel = 1e-10     !relative accuracy
      band_ratio = 2.**(1./3.)
      ! Module AirfoilParams
      Mach = real(U  / p%SpdSound)
  
negative.code:         band_width = freq(i_omega)*(sqrt(band_ratio)-1./sqrt(band_ratio)) * 4. * pi
          ISSUCTION = .TRUE.
          IF (Cf(1) .GT. 0.) THEN
              CALL solve_qk61(omega,f_int2,a,b,answer,abserr,resabs,resasc)   ! wrapper for qk61 routine from slatec library
              Spectrum = D/(4.*pi*R**2.)*answer
              SPL_suction = 10.*log10(Spectrum*DBARH/2.e-5/2.e-5)
              SPLS(i_omega) = SPL_suction + 10.*log10(band_width)
  