Hello,
Triggered by a recent thread ( http://lists.puredata.info/pipermail/pd-dev/2011-07/017022.html) I've written alternative code for Pd classes like phasor~, osc~, and vcf~, in preparation for double precision builds of Pd. These classes currently employ Hoeldrich's method, a fascinating type punning trick for optimization of phase-wrapping jobs. Considering available data types in C, this method can not be translated to double precision input/output format.
Phasor~ and osc~ are typically used by the dozens in Pd setups, and even a modest performance loss could be a show stopper for setups which are on the limit of acceptable CPU load. Fortunately it was possible to define double-ready versions without performance loss, by tedious elimination of redundancy instead of fancy bit hacking. How often must a phase be wrapped? Even at high frequencies like half the Nyquist rate, only one in four iterations goes out of bounds. On average it proves efficient to do the phase wrapping conditionally. Phasor~ profits most, being up to three times faster than before, at moderate frequencies. The others did not speed up so much, but at least none of them is slower than the original version, when both are compared within the same Pd build. Also, the alternative versions do not suffer from phase drift, like the Hoeldrich versions do. I have tested this on OSX / IntelMac for single and double precision builds. Performance may be different on other platforms / architectures. Macro PD_BIGORSMALL was redefined to work with arbitrary precision.
A Pd running with PD_FLOATTYPE double immediately shouted at me that there's a lot more work to do. The audio API (PortAudio in this case) couldn't handle 64 bit output samples from dac~. Vectors are apparently read as type float, and maximum level digital noise is the useless output. I've not delved into this yet. Internally things seems to work well, for what I've seen so far.
Ah, almost forgot to mention the pro's of a double precision Pd, if it will ever work:
- indexing of large audio buffers with ample resolution for interpolation - accurate sine and sinc of small angles, useful for things like fractional delay - FFT with frames larger than 2^17, I hope - long sine sweeps without artifacts - probably many more of such meticulous dsp jobs, you name them
If you're interested in double precision Pd, please find the attached pd_doubletest.zip with C code and test-patches. Let me know if you have test results on different platforms, or if you find stupid bugs in my code. Is anyone working on other aspects of double precision Pd? I'd like to join forces.
Katja Vetter
Hey Katya,
I'm very happy you're working on this, I think its a big and very valuable step for Pd for many reasons. For me, things like accessing large arrays and also working with UNIX timestamps and other large integers directly, make Pd a lot easier in cases that touch on those limitations. It brings Pd 99% to the goal of having a generic "number" type that handles all the floats and ints that are used with any frequency.
Sounds like you've done a fair amount of testing. I think the big question that needs to be answered before this gets included is: can this be done without majoring impacting 32-bit operation? It sounds like you've covered that already. As for 64-bit floats to output, a quick hack to get things working is to just hammer samples down to 32- bits...
.hc
On Jul 27, 2011, at 3:39 AM, katja wrote:
Hello,
Triggered by a recent thread (http://lists.puredata.info/pipermail/pd-dev/2011-07/017022.html ) I've written alternative code for Pd classes like phasor~, osc~, and vcf~, in preparation for double precision builds of Pd. These classes currently employ Hoeldrich's method, a fascinating type punning trick for optimization of phase-wrapping jobs. Considering available data types in C, this method can not be translated to double precision input/output format.
Phasor~ and osc~ are typically used by the dozens in Pd setups, and even a modest performance loss could be a show stopper for setups which are on the limit of acceptable CPU load. Fortunately it was possible to define double-ready versions without performance loss, by tedious elimination of redundancy instead of fancy bit hacking. How often must a phase be wrapped? Even at high frequencies like half the Nyquist rate, only one in four iterations goes out of bounds. On average it proves efficient to do the phase wrapping conditionally. Phasor~ profits most, being up to three times faster than before, at moderate frequencies. The others did not speed up so much, but at least none of them is slower than the original version, when both are compared within the same Pd build. Also, the alternative versions do not suffer from phase drift, like the Hoeldrich versions do. I have tested this on OSX / IntelMac for single and double precision builds. Performance may be different on other platforms / architectures. Macro PD_BIGORSMALL was redefined to work with arbitrary precision.
A Pd running with PD_FLOATTYPE double immediately shouted at me that there's a lot more work to do. The audio API (PortAudio in this case) couldn't handle 64 bit output samples from dac~. Vectors are apparently read as type float, and maximum level digital noise is the useless output. I've not delved into this yet. Internally things seems to work well, for what I've seen so far.
Ah, almost forgot to mention the pro's of a double precision Pd, if it will ever work:
- indexing of large audio buffers with ample resolution for
interpolation
- accurate sine and sinc of small angles, useful for things like
fractional delay
- FFT with frames larger than 2^17, I hope
- long sine sweeps without artifacts
- probably many more of such meticulous dsp jobs, you name them
If you're interested in double precision Pd, please find the attached pd_doubletest.zip with C code and test-patches. Let me know if you have test results on different platforms, or if you find stupid bugs in my code. Is anyone working on other aspects of double precision Pd? I'd like to join forces.
Katja Vetter <pd_doubletest.zip>_______________________________________________ Pd-dev mailing list Pd-dev@iem.at http://lists.puredata.info/listinfo/pd-dev
----------------------------------------------------------------------------
"It is convenient to imagine a power beyond us because that means we don't have to examine our own lives.", from "The Idols of Environmentalism", by Curtis White
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1
On 07/27/2011 11:11 PM, Hans-Christoph Steiner wrote:
like you've covered that already. As for 64-bit floats to output, a quick hack to get things working is to just hammer samples down to 32-bits...
i don't think that's such a great idea. loads of problems (mainly with granular synthesis or other applications where you want to access large tables sample accurately in the signal(!) flow) can simply be fixed by making signals be 64bit too.
and then, quite some infrastructure code makes no clear separations between t_float and t_sample, so it might be simpler to make Pd use doubles throughout and not just for one type of numbers.
gfasdrm IOhannes
On Jul 27, 2011, at 5:47 PM, IOhannes m zmölnig wrote:
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1
On 07/27/2011 11:11 PM, Hans-Christoph Steiner wrote:
like you've covered that already. As for 64-bit floats to output, a quick hack to get things working is to just hammer samples down to 32-bits...
i don't think that's such a great idea. loads of problems (mainly with granular synthesis or other applications where you want to access large tables sample accurately in the signal(!) flow) can simply be fixed by making signals be 64bit too.
and then, quite some infrastructure code makes no clear separations between t_float and t_sample, so it might be simpler to make Pd use doubles throughout and not just for one type of numbers.
I'm saying only as the final output stage to portaudio as a temporary hack to get things working. Its not a good idea otherwise.
.hc
----------------------------------------------------------------------------
“We must become the change we want to see. - Mahatma Gandhi
Hi all --
I'm not sure it's done right, but my intention in s_audio_pa.c is to use 'float' when talking to the portaudio API and t_sample to talk to the rest of Pd -- so t_sample could be made double without affecting portaudio.
The only situation I can imagine in which t_sample might want to differ from t_float is to do ficed-point audio... but I think nowadays that's almost never needed.
cheers Miller
On Wed, Jul 27, 2011 at 05:56:01PM -0400, Hans-Christoph Steiner wrote:
On Jul 27, 2011, at 5:47 PM, IOhannes m zmölnig wrote:
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1
On 07/27/2011 11:11 PM, Hans-Christoph Steiner wrote:
like you've covered that already. As for 64-bit floats to output, a quick hack to get things working is to just hammer samples down to 32-bits...
i don't think that's such a great idea. loads of problems (mainly with granular synthesis or other applications where you want to access large tables sample accurately in the signal(!) flow) can simply be fixed by making signals be 64bit too.
and then, quite some infrastructure code makes no clear separations between t_float and t_sample, so it might be simpler to make Pd use doubles throughout and not just for one type of numbers.
I'm saying only as the final output stage to portaudio as a temporary hack to get things working. Its not a good idea otherwise.
.hc
“We must become the change we want to see. - Mahatma Gandhi
Pd-dev mailing list Pd-dev@iem.at http://lists.puredata.info/listinfo/pd-dev
On Wed, Jul 27, 2011 at 11:11 PM, Hans-Christoph Steiner hans@at.or.atwrote:
I think the big question that needs to be answered before this gets
included is:
can this be done without majoring impacting 32-bit operation?
The intention is: - Pd source code with unaltered functionality - compilable with pd floattype 'float' or 'double' - as little conditional compilation as possible - no performance loss respective to current Pd
Based on test results I think it's possible to rewrite the code in such a way that single precision Pd will not be affected in any way. I still have to rewrite tabosc~ which also uses Hoeldrich's method, this will be easy.
As for 64-bit floats to output, a quick hack to get things working is to just hammer samples down to 32-bits...
I was looking for a suitable spot in the code to do this. First looked at dac~, but since there may be many dac~s instantiated this is not most efficient. Then I found sys_send_dacs(), where the integrated sample values are checked for max absolute value. It is however not possible to do a simple typecast here because samples are just stored back into *sys_soundin and *sys_soundout which are type t_sample. Maybe dac~ should integrate samplevalues in an intermediate vector of type t_sample. And then, in sys_send_dacs(), integrated samples could be checked, cast to float and stored into *sys_soundout vector of type float. And something similar for the input. That's what I'll try.
Katja
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1
On 07/28/2011 02:15 AM, katja wrote:
I was looking for a suitable spot in the code to do this. First looked at dac~, but since there may be many dac~s instantiated this is not most efficient. Then I found sys_send_dacs(), where the integrated sample values are checked for max absolute value. It is however not possible to do a simple typecast here because samples are just stored back into *sys_soundin and *sys_soundout which are type t_sample. Maybe dac~ should integrate samplevalues in an intermediate vector of type t_sample. And then, in sys_send_dacs(), integrated samples could be checked, cast to float and stored into *sys_soundout vector of type float. And something similar for the input. That's what I'll try.
imho, the only sensible place to do this is the actual audio backend code (s_audio_*.c), since your audio backend might support double floats (e.g. jackd is prepared for that already, though i don't known whether anybody ever used it that way).
you don't want to do the calculations in double, then convert to single only to output as double, do you?
fgmadr IOhannes
2011/7/28 IOhannes m zmölnig zmoelnig@iem.at wrote:
imho, the only sensible place to do this is the actual audio backend code (s_audio_*.c)
right
your audio backend might support double floats (e.g. jackd is prepared for that already, though i don't known whether anybody ever used it that way).
Thanks for pointing to that. Anyway some sort of galvanic separation between arrays with t_sample and arrays passed from/to audio API must be introduced, in case they need be of different (primitive) type. PortAudio does not currently support 64 bit samples, for one thing. My first concern is to hear decent sound coming out of double precision Pd, that would really make me happy for a while.....
Katja
Managed to get a double precision Pd working with PortAudio, after small modifications in s_audio_pa.c.
In function pa_open_audio(), *soundin and *soundout (pointers to type t_sample) used to be aliased to *pa_soundin resp. *pa_soundout, pointers to type float. If instead, *pa_soundin and *pa_soundout (and it's own aliases) are declared of type t_sample, the hard copying of samples from PortAudion arrays to Pd arrays and vice versa does the typecast between types of different sizes correctly (tested for Pd 0.42.6, I'll do it for 0.43 soon as possible).
Since I am not familiar with this part of Pd code, my proceedings are slow and I've not found the time to scan s_audio_jack.c, s_audio_alsa.c etcetera for similar issues. Yet, on the dsp side of things, I feel we're close to a properly functioning double precision Pd.
In the meantime I'd like to toss another topic of concern. While testing, I found that a single precision Pd can easily load double precision executables if they are present in the search path for some reason. The other way round will undoubtedly be possible as well - the C symbols in single and double precision builds are identical. Pd crashes on such a size mismatch. Some protection against the loading of 'type aliens' must be implemented before double precision Pd (and externals in particular) can be distributed as stable release, whenever that will happen. And this protection must preferably be backwards compatible, in such a way that old Pd installations from the 'single only era' can not accidentally load new double precision classes.
Katja
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1
On 07/30/2011 01:04 PM, katja wrote:
Since I am not familiar with this part of Pd code, my proceedings are slow and I've not found the time to scan s_audio_jack.c, s_audio_alsa.c etcetera for similar issues. Yet, on the dsp side of things, I feel we're close to a properly functioning double precision Pd.
i had correct output using double precision samples using (iirc) ALSA and jack, back then when i started the patch set (which eventually made it into Pd-0.42), so that shouldn't be a problem.
looking at the code, the relevant changes are definitely there for jack; and alsa will most likely work if it's not trying to do mmap() - if the device does support mmap we might have problems.
distributed as stable release, whenever that will happen. And this protection must preferably be backwards compatible, in such a way that old Pd installations from the 'single only era' can not accidentally load new double precision classes.
totally true.
fmgsadr IOhannes
2011/7/30 IOhannes m zmölnig zmoelnig@iem.at wrote:
looking at the code, the relevant changes are definitely there for jack; and alsa will most likely work if it's not trying to do mmap() - if the device does support mmap we might have problems.
In s_audio_alsamm.c there is an implicit typecast from t_sample to float in line 1208 where the samples are copied (function alsamm_send_dacs()). And for the input a typecast is at line 1316. This makes me think doubles could work with alsa mmap as well. Can't test it though.
I will not have time to work on the Pd double precision thing till after Pd con, but will resume testing later, on OSX, Linux and maybe Windows.
Katja
-
Hans-Christoph Steiner -
IOhannes m zmölnig -
katja -
Miller Puckette