anchor.diff:          call log_note('bolometric correction tables successfully initialized')
       end subroutine bc_init
   
  -    subroutine bc_eval(teff, logg, logl, magnitudes)
  +    subroutine bc_eval(teff, logg, logl, mag_length, magnitudes)
           use interpolate, only: lagrange
   
  -        integer :: i, j, k, teff_index, logg_index
  +        integer :: i, j, k, teff_index, logg_index, mag_length
           real(dp) :: m_bol
           real(dp), intent(in) :: teff, logg, logl
           real(dp), parameter  :: m_bol_sun = 4.74
           real(dp), dimension(4) :: coeffs
           real(dp), dimension(:,:), allocatable :: bc_loggs
  -        real(dp), dimension(:), allocatable, intent(out) :: magnitudes
  -
  -!f2py intent(in) teff, logg, logl
  -!f2py intent(out) magnitudes
  -
  -        if (allocated(magnitudes) .eqv. .false.) then
  -            call log_note('allocating memory for magnitudes in bc_eval')
  -            allocate(magnitudes(n_passbands))
  -        else
  -            call log_warn('magnitudes in bc_eval already allocated')
  -            deallocate(magnitudes)
  -            allocate(magnitudes(n_passbands))
  -        end if
  +        real(dp), dimension(mag_length), intent(out) :: magnitudes
  +
  +        !f2py intent(in) teff, logg, logl, mag_length
  +        !f2py intent(out) magnitudes
  +        !f2py depend(mag_length) magnitudes
  +
  +        ! allocatable output arrays don't play well with f2py
  +        !
  +        !if (allocated(magnitudes) .eqv. .false.) then
  +        !    call log_note('allocating memory for magnitudes in bc_eval')
  +        !    allocate(magnitudes(n_passbands))
  +        !else
  +        !    call log_warn('magnitudes in bc_eval already allocated')
  +        !    deallocate(magnitudes)
  +        !    allocate(magnitudes(n_passbands))
  +        !end if
   
           ! define bolometric magnitude
           m_bol = m_bol_sun - 2.5*logl
  
positive.code: 
          real(dp), intent(in) :: feh, afe
  
          character(len=*), intent(in) :: brand
          character(len=*), dimension(:), intent(in) :: filters
  
  !f2py intent(in) feh, afe, brand, filters
  
          n_passbands = size(filters, 1)
          if (allocated(passbands) .eqv. .false.) allocate(passbands(n_passbands))
  
  
negative.code:             passbands(i) = trim(filters(i))
          end do
  
          call log_note('initializing bolometric correction tables')
          ! initialize bolometric correction tables
          select case (trim(bc_type))
              case ('marcs')
  