How to implement a filter associated to a specific wavelet
You've mentioned Butterworth filters for doing the wavelet analysis using bior6.8. If you want to perform the Discrete Wavelet Transform using some specific wavelet, then you must use its Perfect Reconstruction Filter Bank. Each wavelet function has its associated set of filter values for decomposition and reconstruction - they are calculated from the Mother and Daughter wavelet.

Since you were using Python, I am going to give you an example using PyWavelets. It is as simple as creating a wavelet object and asking for the filter bank. Here is a plot of the filter values:



Since these are FIR filter coefficients, then we can calculate the amplitude response by taking the DFT of the coefficient values:





Values of filter coefficients:

dec_lo = [0.0,
 0.0019088317364812906,
 -0.0019142861290887667,
 -0.016990639867602342,
 0.01193456527972926,
 0.04973290349094079,
 -0.07726317316720414,
 -0.09405920349573646,
 0.4207962846098268,
 0.8259229974584023,
 0.4207962846098268,
 -0.09405920349573646,
 -0.07726317316720414,
 0.04973290349094079,
 0.01193456527972926,
 -0.016990639867602342,
 -0.0019142861290887667,
 0.0019088317364812906]

dec_hi = [-0.0,
 0.0,
 -0.0,
 0.014426282505624435,
 -0.014467504896790148,
 -0.07872200106262882,
 0.04036797903033992,
 0.41784910915027457,
 -0.7589077294536541,
 0.41784910915027457,
 0.04036797903033992,
 -0.07872200106262882,
 -0.014467504896790148,
 0.014426282505624435,
 -0.0,
 0.0,
 -0.0,
 0.0]

rec_lo = [0.0,
 0.0,
 0.0,
 0.014426282505624435,
 0.014467504896790148,
 -0.07872200106262882,
 -0.04036797903033992,
 0.41784910915027457,
 0.7589077294536541,
 0.41784910915027457,
 -0.04036797903033992,
 -0.07872200106262882,
 0.014467504896790148,
 0.014426282505624435,
 0.0,
 0.0,
 0.0,
 0.0]

rec_hi = [0.0,
 -0.0019088317364812906,
 -0.0019142861290887667,
 0.016990639867602342,
 0.01193456527972926,
 -0.04973290349094079,
 -0.07726317316720414,
 0.09405920349573646,
 0.4207962846098268,
 -0.8259229974584023,
 0.4207962846098268,
 0.09405920349573646,
 -0.07726317316720414,
 -0.04973290349094079,
 0.01193456527972926,
 0.016990639867602342,
 -0.0019142861290887667,
 -0.0019088317364812906]




A complementary code which might be useful for you to investigate other wavelets:

import pywt
import numpy as np
import matplotlib.pyplot as plt

# Create wavelet and extract the filters
wavelet_name = 'bior6.8'
wavelet = pywt.Wavelet(wavelet_name)
dec_lo, dec_hi, rec_lo, rec_hi = wavelet.filter_bank

# Filter coefficients
plt.figure()
plt.subplot(221)
plt.stem(dec_lo)
plt.grid()
plt.title('{} low-pass decomposition filter'.format(wavelet_name))
plt.subplot(222)
plt.stem(dec_hi)
plt.grid()
plt.title('{} high-pass decomposition filter'.format(wavelet_name))
plt.subplot(223)
plt.stem(rec_lo)
plt.grid()
plt.title('{} low-pass reconstruction filter'.format(wavelet_name))
plt.subplot(224)
plt.stem(rec_hi)
plt.grid()
plt.title('{} high-pass reconstruction filter'.format(wavelet_name))

# Frequency responses
dec_lo_fr = np.abs(np.fft.rfft(dec_lo, 128))
dec_hi_fr = np.abs(np.fft.rfft(dec_hi, 128))
rec_lo_fr = np.abs(np.fft.rfft(rec_lo, 128))
rec_hi_fr = np.abs(np.fft.rfft(rec_hi, 128))

plt.figure()
plt.subplot(211)
plt.plot(dec_lo_fr, label='Low-pass')
plt.hold(True)
plt.plot(dec_hi_fr, label='High-pass')
plt.grid()
plt.legend()
plt.title('Frequency responses of {} decomposition filters'.format(wavelet_name))
plt.subplot(212)
plt.plot(rec_lo_fr, label='Low-pass')
plt.hold(True)
plt.plot(rec_hi_fr, label='High-pass')
plt.grid()
plt.legend()
plt.title('Frequency responses of {} reconstruction filters'.format(wavelet_name))

plt.show()




Regarding the actual filtering process you can either use scipy.signal.lfilter or np.convolve (x being the input signal):

y = lfilter(dec_lo, 1, x)


or 

y = convolve(x, dec_lo)


I would suggest using PyWavelets to do all of that for you, instead of re-inventing the wheel. However, If you really want to implement a single level of DWT on your own then I second the numpy option. Here is a complete example:

t = np.linspace(0, 1.0, 128)
x = np.sin(2*np.pi*10*t)

# Perform manual DWT and decimate
cA = np.convolve(x, dec_lo)[1::2]
cD = np.convolve(x, dec_hi)[1::2]

plt.figure()
plt.subplot(211)
plt.plot(cA)
plt.grid()
plt.title('Approximation coefficients')

plt.subplot(212)
plt.plot(cD)
plt.grid()
plt.title('Detail coefficients')

plt.show()


Which gives: