Image stabilization using OpenCV

This application deals with video of neural recordings. In such recordings, feature sizes are small. On top of it, recordings are quite noisy. Animal head movements introduces sharp shakes. Out of the box video stabilizer may not work very well on such recordings. Though there are quite a lot of plugins for ImageJ to do such a work, I haven’t compared their performance with this application. This application is hosted here and a demo video is available on youtube here .

The summary of basic principle is following:

0. Collect all frames in a list/vector.

1. Use bilateral filter to smooth out each frame. Bilateral filter smoothens image without distorting the edges (well to a certain extent).

2.  Calculate optical flow between previous frame and current frame. This is a proxy for movement. Construct a transformation and store them in a vector. OpenCV function `goodFeatureToTrack` does almost all the work for us.

  1. Take average of these transformations and apply it on each frame of original recording; that’s correct motion.
Posted in Application, Utility | Tagged , , | Leave a comment

A csv reader based on Haskell-cassava library : performance

I implemented my own csv reader using cassava library. The reader from missingh library was taking too long (~ 17 seconds) for a file with 43200 lines. I compared the result with python-numpy and python-pandas csv reader. Below is rough comparison.

cassava (ignore #) 3.3 sec
cassava (no support for ignoring #) 2.7 sec
numpy loadtxt > 10 sec
pandas read_csv 1.5 sec

As obvious, pandas does really well at reading csv file. I was hoping that my csv reader would do better but it didn’t. But it still beats the parsec based reader hands down.

The code is here

Posted in Haskell, Uncategorized | Tagged , , | Leave a comment

Thresholding numpy array, well sort of

Here is a test case

>>> import numpy as np
>>> a = np.array( [ 0.0, 1, 2, 0.2, 0.0, 0.0, 2, 3] )

I want to turn all non-zero elements of this array to 1. I can do it using np.where and numpy indexing.

>>>  a[ np.where( a != 0 ) ] = 1
>>> a
array([ 0.,  1.,  1.,  1.,  0.,  0.,  1.,  1.])

note: np.where returns the indices where the condition is true. e.g. if you want to change all 0 to -1.

>>> a[ np.where[ a == 0] ] = -1.0

That’s it. Checkout np.clip as well.

Posted in Python, Uncategorized | Tagged , | Leave a comment

C++11 library to write std::vector to numpy `npy` file (format 2)

Here is the use scenario: You want to write/append your STL vector to a numpy binary file which you can process using python-numpy/scipy.

Another scenario: Your application generates a huge amount of data. If it is not written to disk, you will end up consuming so much RAM that kernel will kill the application. Occasionally you have to flush the data to a file to keep going.


The quick approach is to write to a text file (which proper formatting of double datatype  — check boost::lexical_cast<double>(string) convenient if already using boost, or sprintf which is fastest, stringstreams are comparatively slow), also std::to_string does not have very high precision). Unfortunately text-files are bulky; and implementing binary format is too much work.

Numpy format is quite minimal It has two version. There is already a pretty good library which supports version 1 of numpy format And it can also write to `npz` format (zipped). The size of headers are limited in version 1.


I needed the ability to write bigger header (numpy structured arrays). I wrote this . It writes/appends data to a given file in numpy file format version 2. See the file for more details. It only supports writing to npy format. That’s all! Unfortunately, only numpy >= 1.9  supports version 2. So to make this library a bit more useful, version 1 is also planned.

Posted in Library | Tagged | Leave a comment

Safe level of mercury in soil, water and food

Following are verbatim from book chapter ( ) .


The World Health Organization guideline value for inorganic mercury vapor is 1 μg/m3 as an annual average.227 A tolerable concentration is 0.2 μg/m3 for long-term inhalation exposure to elemental mercury vapor, and a tolerable intake of total mercury is 2 μg/kg body weight per day.70


Food and Agriculture Organization/World Health Organization Codex Alimentarius—Commission guideline levels for methylmercury in fish 0.5 mg/kg for predatory fish (such as shark, swordfish, tuna, pike and others) is 1 mg/kg.228


For methylmercury, the Joint Food and Agriculture Organization/World Health Organization Expert Committee on Food Additives (JECFA) set in 2004 a tolerable weekly intake of 1.6 μg/kg body weight per week to protect the developing fetus from neurotoxic effects.229 JEFCA230 confirmed this provisional tolerable weekly intake level, taking into account that adults might not be as sensitive as the developing fetus, in 2003 (JECFA/61/SC and 2006 (JECFA/67/SC,232


United Nations Environment Programme Global Mercury Assessment quotes for soil, preliminary critical limits to prevent ecological effects due to mercury in organic soils with 0.07-0.3 mg/kg for the total mercury content in soil.4


Posted in Uncategorized | Tagged | Leave a comment

Benchmark ODE solver: GSL V/s Boost Odeint library

For our neural simulator, MOOSE, we use GNU Scientific Library (GSL) for random number generation, for solving system of non-linear equations, and for solving ODE system.

Recently I checked the performance of GSL ode solver V/s Boost ode solver; both using Runge-Kutta 4. The boost-odeint outperformed GSL by approximately by a factor of 4. The numerical results were same. Both implementation were compiled with -O3 switch.

Below are the numerical results. In second subplot, a molecule is changing is concentration ( calcium ) and on the top, another molecule is produced (CaMKII). This network has more than 30 reactions and approximately 20 molecules.

GSL took approximately 39 seconds to simulate the system for 1 year. While Boost-odeint took only 8.6 seconds!


Posted in Biological systems, Numerical computation, Uncategorized | Tagged , , , , | Leave a comment

Multi-dimentional root fiding using Newton-Raphson method

In short, It does the same thing as GNU-gsl `gnewton` solver for finding roots of a system of non-linear system.

The standard caveats apply when finding the multi-dimensional roots.

The templated class (must be compiled with  -std=c++11). It uses boost::numeric::ublas library.

The code is available here:

The header file contains the class, file `test_root_finding.cpp` shows how to use it.


Posted in Uncategorized | Tagged , , , , | Leave a comment