Talk:Discrete sine transform

"(since the Fourier transform of a real and odd function is real and imaginary)"

Shouldn't that be "imaginary and odd"?

Right, whoops. —Steven G. Johnson 19:33, 20 December 2005 (UTC)[reply]

Computing DST using DCT

I tested this in Matlab and came to the conclusion, that the article is wrong regarding this topic. According to my test the order doesn`t matter for DST-IV. "signflip - DCT-IV - reverse" works as well as "reverse - DCT-IV - signflip". Furthermore, for DST-II / DST-III it should be the other way around than described in the article, i.e for DST-II: "signflip - DCT-II - reverse" and for DST-III: "reverse - DCT-III - signflip".

Operating on real input

Does anybody know why the DST (and the DCT) is the defined as only operating on real numbers? I realize that it's more often used with real than with complex data, but I can see no particular reason why it should not apply to complex numbers as well, with complex output -- in fact, I just needed a complex, fast DST, and almost all references I've found dealt with real DSTs only.

Of course, any linear operation on real numbers also trivially operates on complex numbers. The point here is that the transform itself corresponds to a purely real matrix. —Steven G. Johnson

Apart from that, I really miss the formulas for explicit factorization and conversion to a DFT of the same size; after much pain, I've managed to figure out the ones for DST-II (factorization and DFT of the same size (complex -> complex, missing real -> real but it is much easier to find in the literature) and DST-III (DFT of the same size only, complex -> complex); I've been unable to find a good form to present the material in encyclopedic form, though, but if anybody want a pre-copy of the relevant paper (it's for a semester project, so it's not done/public yet), feel free to contact me and I'll send you the calculations. - Sesse 21:09, 16 October 2006 (UTC)[reply]

The only good algorithms for expressing a DFT in terms of a (real-input) DFT of the same size are for the DST-II/III, as described by Makhoul in 1980. For the DST-I there is an algorithm described in e.g. Numerical Recipes, but it has accuracy problems (errors growing as √N instead of √log N). For the DST-IV there was an algorithm described by Chan and Ho in 1990, but it has similar accuracy problems. To perform the DST-I and DST-IV accurately, as far as I know, requires more complicated algorithms that involve e.g. breaking them down into transforms of half the size. See e.g. www.fftw.org for source code and references. —Steven G. Johnson 01:25, 17 October 2006 (UTC)[reply]

Anyway, the algorithms are probably too complicated to derive on this page; at some point they might get their own articles, I suppose. One point that the article might want to make is that it is trivial to get an algorithm for the DST II-IV from an algorithm for the DCT II-IV, just by reversing the input and flipping the sign of every other output (or vice versa for type III). I learned this trick from FFTPACK; I'm not sure how widely known it is but it's quite useful since DCT algorithms are much more common than DST algorithms. —Steven G. Johnson 01:31, 17 October 2006 (UTC)[reply]

Well, they could at least be described; a full proof isn't neccessary, but after looking at papers for twelve hours straight you start to tire a bit of the "it's possible, but we don't really care to tell you exactly how" mantra :-) I actually looked at Makhoul's article and tried to reproduce the calculations for the FST (in vain), eventually finding what I was looking for in the FFTW3 source code (I guess you wrote those comments -- thanks! :-) I had almost given up...), and was able to figure it out for the DST. (I was a bit surprised to see the FFTw3 source code refer to Makhoul's article, though; the extension explanation given there didn't match my understanding of the article at all.) I didn't use the sign flipping trick, though, I used a variant with a complex multiplication instead; it simply wasn't a very relevant optimization for my use case. -Sesse 00:49, 19 October 2006 (UTC)[reply]

I agree that the Makhoul article didn't derive it via extension into a DFT of length 4N; I actually only came across the Makhoul paper long after deriving the algorithm for myself with the help of the FFTPACK source code. Unfortunately, if you look at the literature, there are countless papers deriving DCT algorithms by painful manual trigonometry...and if I look closely at them I've invariably found that they are simply a standard complex FFT algorithm (usually Cooley-Tukey) with the redundant operations thrown out for the appropriate real-symmetric DFT. Someday, I would like to convince people in the field that it's crazy to re-derive the same algorithm 16 times. But in the meantime, it's important for scholarly and historical purposes to cite the paper that first figured out an algorithm, even if they derived it in a way you don't care for. —Steven G. Johnson 06:51, 19 October 2006 (UTC)[reply]

Yes, I've looked at the papers, and felt the pain :-) Anyhow, that clears it up for me; I'll be sure to get my references right. BibTeX doesn't really have a citation type for "comments in source code", though :-) -Sesse 10:50, 19 October 2006 (UTC)[reply]

Request check for DST-I formula

I can't verify this myself, but my intuition tells me that because the indices are ranging from 0 to N-1, then the part with π/N+1 should instead have π/N. Compare with Discrete Fourier transform#Definition. If not, I would like to know why, as it would help in learning this bit of math. LokiClock (talk) 13:57, 22 December 2009 (UTC)[reply]

Nope, it's correct as-is. The basic reason is that a DST-I of, say, 3 data points (a,b,c) is equivalent to a DFT of 8 data points (0,a,b,c,0,-c,-b,-a), a real sequence with odd symmetry. In the DFT, 2π/8 appears, so in the DST-I a π/4 = π/(N+1) appears. I added this clarification to the article. — Steven G. Johnson (talk) 19:19, 18 February 2010 (UTC)[reply]

That makes sense, thanks. LokiClock (talk) 14:27, 19 February 2010 (UTC)[reply]