Ron Kneusel emailed our internal PWUISATAI (People Who Use IDL Sitting Around Talking About IDL) group yesterday with a clever function he saw in some CUDA C code to calculate the next power of two greater than or equal to a given integer. His IDL translation:
function next_power_2, x compile_opt idl2, logical_predicate n = x - 1 ; protects input n = ishft(n, -1) or n n = ishft(n, -2) or n n = ishft(n, -4) or n n = ishft(n, -8) or n n = ishft(n,-16) or n n = ishft(n,-32) or n n = ishft(n,-64) or n return, ++n end
Both ISHFT and the OR operator are used here to perform bitwise operations on the input integer.
Here’s an example of using NEXT_POWER_2:
IDL> a = 3565946L IDL> b = next_power_2(a) IDL> print, b 4194304 IDL> factor, b 22 4194304 = 2
I’ve used Ray Sterner’s FACTOR, a part of the IDL astrolib, to learn that 4194304 is 222.
Atle Borsholm replied to the group with an alternate take:
function n2, x return, ishft(1ull,total(ishft(1ull,indgen(64)) lt x, /integer)) end
that gives the same result:
IDL> b = n2(a) IDL> print, b 4194304
Atle commented that although his version is shorter, Ron’s may be faster. I was curious, so I made a simple time test:
pro time_test_next_power_2 compile_opt idl2 n_iter = 1e6 x = 275259L t0 = systime(/seconds) for i=1, n_iter do !null = next_power_2(x) t1 = systime(/seconds) print, 'NEXT_POWER_2:', t1-t0, format='(a15,f12.8,1x,"s")' t0 = systime(/seconds) for i=1, n_iter do !null = n2(x) t1 = systime(/seconds) print, 'N2:', t1-t0, format='(a15,f12.8,1x,"s")' end
and ran it on my laptop:
IDL> time_test_next_power_2
NEXT_POWER_2: 2.14159489 s
N2: 2.77271295 s
So, Ron’s version is slightly faster. In either case, we hope you might find these functions handy!
Update: An even simpler version from a zero-padding routine I’d written long ago:
IDL> a = 3565946L IDL> b = 2UL^ceil(alog(a)/alog(2)) IDL> print, b 4194304
Note: I’ll be out on vacation for a bit, so I have some guest bloggers lined up for the next few weeks.