Mercurial > sdl-ios-xcode
comparison src/audio/SDL_audiocvt.c @ 2659:8da698bc1205 gsoc2008_audio_resampling
Fixed lots of bugs in FIR filtering. Fixed point code is closer to working, but there seems to be overflow in the FIR filter resulting in distortion.
author | Aaron Wishnick <schnarf@gmail.com> |
---|---|
date | Sun, 22 Jun 2008 00:36:35 +0000 |
parents | de29a03cb108 |
children | a55543cef395 |
comparison
equal
deleted
inserted
replaced
2658:de29a03cb108 | 2659:8da698bc1205 |
---|---|
52 #endif | 52 #endif |
53 /* Confirmed that SDL_FixMpy16 works, need to check 8 and 32 */ | 53 /* Confirmed that SDL_FixMpy16 works, need to check 8 and 32 */ |
54 #define SDL_FixMpy16(a, b) ((((Sint32)a * (Sint32)b) >> 14) & 0xffff) | 54 #define SDL_FixMpy16(a, b) ((((Sint32)a * (Sint32)b) >> 14) & 0xffff) |
55 #define SDL_FixMpy8(a, b) ((((Sint16)a * (Sint16)b) >> 7) & 0xff) | 55 #define SDL_FixMpy8(a, b) ((((Sint16)a * (Sint16)b) >> 7) & 0xff) |
56 /* Everything is signed! */ | 56 /* Everything is signed! */ |
57 #define SDL_Make_1_7(a) (Sint8)(a * 128.0f) | 57 #define SDL_Make_1_7(a) (Sint8)(a * 127.0f) |
58 #define SDL_Make_1_15(a) (Sint16)(a * 32768.0f) | 58 #define SDL_Make_1_15(a) (Sint16)(a * 32767.0f) |
59 #define SDL_Make_1_31(a) (Sint32)(a * 2147483648.0f) | 59 #define SDL_Make_1_31(a) (Sint32)(a * 2147483647.0f) |
60 #define SDL_Make_2_6(a) (Sint8)(a * 64.0f) | 60 #define SDL_Make_2_6(a) (Sint8)(a * 63.0f) |
61 #define SDL_Make_2_14(a) (Sint16)(a * 16384.0f) | 61 #define SDL_Make_2_14(a) (Sint16)(a * 16383.0f) |
62 #define SDL_Make_2_30(a) (Sint32)(a * 1073741824.0f) | 62 #define SDL_Make_2_30(a) (Sint32)(a * 1073741823.0f) |
63 | 63 |
64 /* Effectively mix right and left channels into a single channel */ | 64 /* Effectively mix right and left channels into a single channel */ |
65 static void SDLCALL | 65 static void SDLCALL |
66 SDL_ConvertMono(SDL_AudioCVT * cvt, SDL_AudioFormat format) | 66 SDL_ConvertMono(SDL_AudioCVT * cvt, SDL_AudioFormat format) |
67 { | 67 { |
1524 #undef iir_fix | 1524 #undef iir_fix |
1525 } | 1525 } |
1526 | 1526 |
1527 /* Apply the windowed sinc FIR filter to the given SDL_AudioCVT struct */ | 1527 /* Apply the windowed sinc FIR filter to the given SDL_AudioCVT struct */ |
1528 static void SDL_FilterFIR(SDL_AudioCVT * cvt, SDL_AudioFormat format) { | 1528 static void SDL_FilterFIR(SDL_AudioCVT * cvt, SDL_AudioFormat format) { |
1529 int n = cvt->len_cvt / (SDL_AUDIO_BITSIZE(format) / 4); | 1529 int n = 8 * cvt->len_cvt / SDL_AUDIO_BITSIZE(format); |
1530 int m = cvt->len_sinc; | 1530 int m = cvt->len_sinc; |
1531 int i, j; | 1531 int i, j; |
1532 | 1532 |
1533 /* Note: this makes use of the symmetry of the sinc filter. | 1533 /* Note: this makes use of the symmetry of the sinc filter. |
1534 We can also make a big optimization here by taking advantage | 1534 We can also make a big optimization here by taking advantage |
1535 of the fact that the signal is zero stuffed, so we can do | 1535 of the fact that the signal is zero stuffed, so we can do |
1536 significantly fewer multiplications and additions. | 1536 significantly fewer multiplications and additions. |
1537 */ | 1537 */ |
1538 #define filter_sinc(type, mult) { \ | 1538 #define filter_sinc(type, mult) { \ |
1539 type *sinc = (type *)cvt->coeff; \ | 1539 type *sinc = (type *)cvt->coeff; \ |
1540 type *state = (type *)cvt->state_buf; \ | 1540 type *state = (type *)cvt->state_buf; \ |
1541 type *buf = (type *)cvt->buf; \ | 1541 type *buf = (type *)cvt->buf; \ |
1542 for(i = 0; i < n; ++i) { \ | 1542 for(i = 0; i < n; ++i) { \ |
1543 if(cvt->state_pos == m) cvt->state_pos = 0; \ | 1543 cvt->state_pos++; \ |
1544 if(cvt->state_pos >= m) cvt->state_pos = 0; \ | |
1545 state[cvt->state_pos] = buf[i]; \ | |
1544 buf[i] = 0; \ | 1546 buf[i] = 0; \ |
1545 for(j = 0; j < m; ++j) { \ | 1547 for(j = 0; j < m; ++j) { \ |
1546 buf[i] += mult(state[j], sinc[j]); \ | 1548 buf[i] += mult(state[(cvt->state_pos - j) % m], sinc[j]); \ |
1547 } \ | 1549 } \ |
1548 } \ | 1550 } \ |
1549 } | 1551 } |
1550 | 1552 |
1551 /* If it's floating point, we don't need to do any shifting */ | 1553 /* If it's floating point, we don't need to do any shifting */ |
1599 if( cvt->coeff == NULL ) { | 1601 if( cvt->coeff == NULL ) { |
1600 return -1; | 1602 return -1; |
1601 } | 1603 } |
1602 | 1604 |
1603 /* Set the length */ | 1605 /* Set the length */ |
1604 cvt->len_sinc = m; | 1606 cvt->len_sinc = m + 1; |
1605 | 1607 |
1606 /* Allocate the floating point windowed sinc. */ | 1608 /* Allocate the floating point windowed sinc. */ |
1607 fSinc = (float *)malloc(m * sizeof(float)); | 1609 fSinc = (float *)malloc(m * sizeof(float)); |
1608 if( fSinc == NULL ) { | 1610 if( fSinc == NULL ) { |
1609 return -1; | 1611 return -1; |
1610 } | 1612 } |
1611 | 1613 |
1612 /* Set up the filter parameters */ | 1614 /* Set up the filter parameters */ |
1613 fc = (cvt->len_mult > cvt->len_div) ? 0.5f / (float)cvt->len_mult : 0.5f / (float)cvt->len_div; | 1615 fc = (cvt->len_mult > cvt->len_div) ? 0.5f / (float)cvt->len_mult : 0.5f / (float)cvt->len_div; |
1616 fc = 0.04f; | |
1614 two_pi_fc = 2.0f * M_PI * fc; | 1617 two_pi_fc = 2.0f * M_PI * fc; |
1615 two_pi_over_m = 2.0f * M_PI / (float)m; | 1618 two_pi_over_m = 2.0f * M_PI / (float)m; |
1616 four_pi_over_m = 2.0f * two_pi_over_m; | 1619 four_pi_over_m = 2.0f * two_pi_over_m; |
1617 m_over_two = (float)m / 2.0f; | 1620 m_over_two = (float)m / 2.0f; |
1618 norm_sum = 0.0f; | 1621 norm_sum = 0.0f; |
1624 fSinc[i] = sinf(two_pi_fc * ((float)i - m_over_two)) / ((float)i - m_over_two); | 1627 fSinc[i] = sinf(two_pi_fc * ((float)i - m_over_two)) / ((float)i - m_over_two); |
1625 /* Apply blackman window */ | 1628 /* Apply blackman window */ |
1626 fSinc[i] *= 0.42f - 0.5f * cosf(two_pi_over_m * (float)i) + 0.08f * cosf(four_pi_over_m * (float)i); | 1629 fSinc[i] *= 0.42f - 0.5f * cosf(two_pi_over_m * (float)i) + 0.08f * cosf(four_pi_over_m * (float)i); |
1627 } | 1630 } |
1628 norm_sum += abs(fSinc[i]); | 1631 norm_sum += abs(fSinc[i]); |
1632 printf("%f\n", fSinc[i]); | |
1629 } | 1633 } |
1630 | 1634 |
1631 /* Now normalize and convert to fixed point. We scale everything to half the precision | 1635 /* Now normalize and convert to fixed point. We scale everything to half the precision |
1632 of whatever datatype we're using, for example, 16 bit data means we use 8 bits */ | 1636 of whatever datatype we're using, for example, 16 bit data means we use 8 bits */ |
1633 | 1637 |
1634 #define convert_fixed(type, fix) { \ | 1638 #define convert_fixed(type, fix) { \ |
1635 norm_fact = 1.0f / norm_sum; \ | 1639 norm_fact = 0.9f / norm_sum; \ |
1640 norm_fact = 0.15f; \ | |
1636 type *dst = (type *)cvt->coeff; \ | 1641 type *dst = (type *)cvt->coeff; \ |
1637 for( i = 0; i <= m; ++i ) { \ | 1642 for( i = 0; i <= m; ++i ) { \ |
1638 dst[i] = fix(fSinc[i] * norm_fact); \ | 1643 dst[i] = fix(fSinc[i] * norm_fact); \ |
1644 printf("%f = 0x%x\n", fSinc[i], dst[i]); \ | |
1639 } \ | 1645 } \ |
1640 } | 1646 } |
1641 | 1647 |
1642 /* If we're using floating point, we only need to normalize */ | 1648 /* If we're using floating point, we only need to normalize */ |
1643 if(SDL_AUDIO_ISFLOAT(format) && SDL_AUDIO_BITSIZE(format) == 32) { | 1649 if(SDL_AUDIO_ISFLOAT(format) && SDL_AUDIO_BITSIZE(format) == 32) { |
1659 break; | 1665 break; |
1660 } | 1666 } |
1661 } | 1667 } |
1662 | 1668 |
1663 /* Initialize the state buffer to all zeroes, and set initial position */ | 1669 /* Initialize the state buffer to all zeroes, and set initial position */ |
1664 memset(cvt->state_buf, 0, cvt->len_sinc * SDL_AUDIO_BITSIZE(format) / 4); | 1670 //memset(cvt->state_buf, 0, cvt->len_sinc * SDL_AUDIO_BITSIZE(format) / 4); |
1665 cvt->state_pos = 0; | 1671 cvt->state_pos = 0; |
1666 | 1672 |
1667 /* Clean up */ | 1673 /* Clean up */ |
1668 #undef convert_fixed | 1674 #undef convert_fixed |
1669 free(fSinc); | 1675 free(fSinc); |
1730 } | 1736 } |
1731 | 1737 |
1732 cvt->len_cvt *= cvt->len_mult; | 1738 cvt->len_cvt *= cvt->len_mult; |
1733 | 1739 |
1734 // Step 2: Use either a windowed sinc FIR filter or IIR lowpass filter to remove all alias frequencies | 1740 // Step 2: Use either a windowed sinc FIR filter or IIR lowpass filter to remove all alias frequencies |
1735 SDL_FilterIIR( cvt, format ); | 1741 SDL_FilterFIR( cvt, format ); |
1736 | 1742 |
1737 // Step 3: Discard unnecessary samples | 1743 // Step 3: Discard unnecessary samples |
1738 #ifdef DEBUG_CONVERT | 1744 #ifdef DEBUG_CONVERT |
1739 printf("Discarding samples by a factor of %u\n", cvt->len_div); | 1745 printf("Discarding samples by a factor of %u\n", cvt->len_div); |
1740 #endif | 1746 #endif |
1853 rate_gcd = SDL_GCD(src_rate, dst_rate); | 1859 rate_gcd = SDL_GCD(src_rate, dst_rate); |
1854 cvt->len_mult = dst_rate / rate_gcd; | 1860 cvt->len_mult = dst_rate / rate_gcd; |
1855 cvt->len_div = src_rate / rate_gcd; | 1861 cvt->len_div = src_rate / rate_gcd; |
1856 cvt->len_ratio = (double)cvt->len_mult / (double)cvt->len_div; | 1862 cvt->len_ratio = (double)cvt->len_mult / (double)cvt->len_div; |
1857 cvt->filters[cvt->filter_index++] = SDL_Resample; | 1863 cvt->filters[cvt->filter_index++] = SDL_Resample; |
1858 SDL_BuildIIRLowpass(cvt, dst_fmt); | 1864 //SDL_BuildIIRLowpass(cvt, dst_fmt); |
1865 SDL_BuildWindowedSinc(cvt, dst_fmt, 12); | |
1859 | 1866 |
1860 /*cvt->rate_incr = 0.0; | 1867 /*cvt->rate_incr = 0.0; |
1861 if ((src_rate / 100) != (dst_rate / 100)) { | 1868 if ((src_rate / 100) != (dst_rate / 100)) { |
1862 Uint32 hi_rate, lo_rate; | 1869 Uint32 hi_rate, lo_rate; |
1863 int len_mult; | 1870 int len_mult; |