forkmap, threadmap

I have an 8 core machine, so vanilla Python increasingly annoys me by running roughly 200*8 times slower than optimized C code (I use OpenMP), rather than 200 times slower than C.

To use more cores, I occasionally use these parallel map functions in Python:

• threadmap
  
• forkmap
  

threadmap spawns a number of threads equal to the number of cores, and distributes work across threads. forkmap requires Cygwin, Mac, or Unix, and does the same thing with fork() and pipe(), basically to get around GIL restrictions in regular Python. Probably a better solution is to use a Python implementation such as IronPython which doesn't have the GIL, but the FFIs for CPython libraries seem spotty in these other Python implementations.

Merely using parallel map() functions is rather limiting; OpenMP-style syntax should probably be added to Python at some point, after the GIL is removed. (As the number of processor cores goes to infinity, either the GIL will be removed or Python will become irrelevant, as languages with similar syntax and fine-grain locking will take over.)

There are libraries which implement parallelism by using spawn() instead of fork() to get Windows compatibility without Cygwin, but I couldn't figure out how exactly these libraries were spawning new copies of my code, so due to lack of documentation on this point I used other solutions.

C++ foreach

I can never remember the C++ for-each loop syntax. Thus, a helper macro I use:

#define for_each(type, it, L) for (type::iterator (it) = (L).begin(); (it) != (L).end(); ++(it))

If L has type T, usage is:

for_each(T, it, L) {
  printf("%d\n", (*it));
}

If you're only using gcc, one can use __typeof__ [1] to simplify this.

RK45 in Python

Runge-Kutta 4th and 5th order adaptive ODE integrator. Generally the scipy integrators will be easier to use, unless you specifically need RK45.
(Code)

Filter numpy images with FFT, Python

Generic linear filter support is not currently built into the Python Imaging Library. This module lets you filter a numpy array against an arbitrary kernel:

>>> I = numpy.asarray(Image.open('test.jpg'))
>>> I = filter(I, [[-1,0,1],[-2,0,2],[-1,0,1]])/8+128
>>> Image.fromarray(numpy.uint8(I)).show()


The Fast Fourier Transform (FFT) is used. The FFT routine included with numpy isn't particularly fast (c.f. FFTW [1]), and in any case using the transform isn't as efficient as applying the filter naively for small filter sizes. Also, for separable kernels (e.g. the Gaussian kernel), it is often faster to perform two 1D convolutions in sequence. However, for large images and filters, using an FFT for convolution is often faster than other approaches. (If an image and filter contain a total of N pixels, then this algorithm takes O(NlogN) time, which is the fastest known time complexity algorithm for the general problem.) Public domain.

Update: A number of image filtering operations are provided in numarray (in the module numpy.numarray.nd_image), such as filtering, image morphology, distance transforms, and segmentation; see [2].

(Source)

framecap - Capture camera input with Allegro, C++

Requires OpenCV, Allegro, and pthreads. Captures frames from a webcam. Usage:

CAPTURE *cap = create_capture();
BITMAP *bmp = capture_frame(cap);

Use capture_frame_nb to capture a frame in nonblocking mode (a NULL pointer is returned if no frame is available). Should work on Mac OS X, Linux, Windows; only tested on Windows.

(Source)

pfm - Read/write 1 and 3 channel PFM files

PFM files are mxnx1 or mxnx3 arrays of floats. C++ code to read/write the format (public domain).
(Source)

OpenCV camera and AVI input

To convert a video to a format that OpenCV can read:

  $ mencoder in.avi -ovc raw -vf format=i420 -o out.avi

This is necessary because OpenCV uses platform specific video libraries, and the intersection of their supported formats is DIB/I420/IYUV. To capture input from a camera or an AVI file, use the following C++ code: [1]. To capture frames in a nonblocking manner in Python, use: [2] with Gary Bishop's cvtypes [3].

Burrows-Wheeler transform in Python

def bwt(s):
  s = s + '\0'
  return ''.join([x[-1] for x in
         sorted([s[i:] + s[:i] for i in range(len(s))])])

def ibwt(s):
  L = [''] * len(s)
  for i in range(len(s)):
    L = sorted([s[i] + L[i] for i in range(len(s))])
  return [x for x in L if x.endswith('\0')][0][:-1]

The Burrows-Wheeler transform [1] is an invertible mathematical transformation which tends to clump related sequences of data together; hence it is useful as a pre-pass for some compression algorithms. The BWT domain in the above code excludes the null character; to include the null character, modify the transformation to also output the row of sorted([s[i:] + s[:i]...]) on which s occurs, then in IBWT simply return this row after the for loop. The Python code above is not particularly fast.

How to put a video on the Web

Install FFMPEG [1], use ffmpeg -i in.avi out.swf, then place <embed src="out.swf" width="400" height="300" loop="false" quality="high"></embed> in your HTML. Optionally, use the -s argument after -i in FFMPEG to change the resolution. If you want playback controls, output to a .flv file instead, download FlowPlayer [1], and modify their example HTML to point to your .flv video.

Ubuntu Linux postinstall notes

After installing an operating system, there are always a zillion changes I have to make to get the operating system operational. I'm new to Linux, so I thought I'd post up a list of changes I made after installing Ubuntu to help out others. These are mostly basic changes such as enabling all software packages, choosing a fast repository, mounting FAT32 drives, making Ctrl+Alt+Del bring up a task manager, installing multimedia software, changing the system DPI, trying to disable the trash feature, and adding an Open Terminal option to the desktop.
(Document)

Eigenvectors of 3x3 symmetric matrix

A C++ source and header file to compute eigenvectors/values of a 3x3 symmetric matrix. Potentially easier than installing EISPACK, LAPACK, or Gandalf if you only need this single function. Takes about 6000 clock cycles per call on my Pentium 4. Public domain. [1]

Incremental flash backup tool

I recently decided that I would back up my documents to a high capacity flash disk. Thus I wrote a backup tool.
Features:

• Dead simple. Something you'd write in an hour, so it won't mess you up by having volumes of code, weird interfaces, complicated man pages, or obscure back up formats.
  
• The flash disk is a mirrored filesystem, so no special tools are required to read the backup.
  
• The existing backed-up filesystem is updated if the backup command is re-issued.
  
• Specified directories are recursively copied; hidden files are not copied. Thus one can flag large files that should not be backed up: make them hidden.
  
• Eventually needs support for checksums, but I'm not too worried at the moment.
  

The program is public domain. Download and run it. Copy its example backup command line into a shell script. In my case, I use something similar to:

  % python backup.py D:\ D:\Documents D:\Code G:\Backup

(Source Code)

Intelligent Minesweeper

A game of Minesweeper, with moves suggested by the entail module I posted earlier. The entail module contains a propositional logic inference algorithm: it lets one determine whether any given expression of binary variables follows logically from a set of assumptions:

>>> s = AssumptionSet(['a => b => c <=> a'])
>>> s.implies('a <=> b')
True
>>> s.assume('c')
>>> s.implies('a')
True

It is rather uncanny how complicated mathematical relationships can be derived ex nihilo using only the entail module. For the Minesweeper game, when the user chooses "Suggest Move," I add only the knowledge that the player possesses to an AssumptionSet, and then suggest a move if the logic algorithm can deduce that any given square is free of a mine. The result is a perfect deterministic minesweeper AI: if any given move on the board can be made safely, then the AI will suggest that square. If no move can be made safely, then no suggestion will be made. Note that this is also an experiment with using HTML widgets in lieu of using a widget toolkit such as wxPython, so bear with the strange appearance of the window.

(Screenshot) (EXE) (Source) (Module entail) (All public domain).

Real-time hair animation

I used the Verlet integration and relaxation method described by Thomas Jakobsen in "Advanced Character Physics" [1] coupled with the animation techniques suggested by Chang, Jin, and Yu in [2] to make a fast hair simulation. My resulting video [3] was captured in real-time. The program can actually be sped up significantly by placing bounding boxes around individual quads, vertices, and segments; however, I ran into bugs when I tried to do this, so I left that optimization out for now. The results are interesting because they indicate that in video games, hair can be animated in real-time with relative ease. The source code has been placed in the public domain.
(Video) (Writeup) (Source Code)

Automatic Python imports with autoimp

I got sick of writing "import X" in Python. So, I created the public domain module autoimp, which imports all modules automatically:

>>> from autoimp import *
>>> os.stat('.')
>>> Image.open('test.bmp')
>>> pylab.plot([1,2],[3,4])
>>> scipy.linalg.eig([[1,2],[3,4]])
>>> ...

Thus one no longer needs to write "import X". It would take too long to load every module when one writes "from autoimp import *", so the imported modules are actually proxy objects which lazily load when they are first used.
(Module Website)

Programming Language Productivity

Norvig [1] commented on two productivity studies by Prechelt [2] and Garret [3]. The results for the productivity of the languages isn't immediately apparent, so I extracted the relevant information as a "productivity summary" [4]. The results strongly indicate that programmers who solve a problem with fewer lines of code are more productive, irrespective of the language used. Thus languages which are more succinct could be said to be more productive.

For entertainment purposes, I solved the given problem in Python; it took me 1 hour 10 minutes, with 31 lines of code. I spent 10 minutes reading the problem, and 40 minutes fixing mistakes that I had made due to not reading the problem carefully enough. There is a lesson here...

(Productivity Summary) (My Solution)

Eight queens in Python

An elegant way to find all solutions to the 8-queens problem: generate all 8-castles solutions by a permutations() generator and filter the diagonals by checking that the set of diagonals represented has 8 elements. I implemented this algorithm in Python (17 lines) and C++ (66 lines).
(Python Code) (C++ Code)

Shedskin benchmarks

Shedskin [1] is an experimental Python to C++ compiler. I benchmarked it for some of the microbenchmarks in The Computer Language Shootout [2]. The results are interesting: for some benchmarks, the compiled Python code runs at nearly the speed of C.
(These benchmarks are out of date, so I removed the link).

rencode -- Reduced length encodings

A space-efficient Python serialization module based on bencode,
which can be used to decode strings from untrusted sources. Works well for complex, heterogeneous data structures with many small elements. The encodings take considerably less space than bencode, gherkin, and Recipe 415503 [1] for that use case, and this module is faster, too. (Source code), (Comparison).

Allegro Quick Reference in Japanese

The Allegro Quick Reference guide that I wrote in 2002 [1] has been translated into Japanese [2] by Ohsawa Hirotaka. Cool.