Type-safe scripting with C++ (and other weird explorations)

The why and the what(-the-f***)

Let’s start this short tale with some background. 

For reasons unclear, I’ve started working on my n-th abandonable side-project. Much detail isn’t necessary here: it’s basically a C++ library for performing simple math operations— averaging, sum, standard deviation, autocorrelation and the like. Each operation is implemented as Functor: instantiate it, call it’s operator(), and done! For a concrete example, consider the SumFunctor below (ignore the inheritance for now), which offers overloads to sum two integers or sum (concatenate) two strings

However, there’s a catch: I also want to be able to invoke these functors from some sort of scripting language. The exact language isn’t relevant (though I have my eye on Tcl); suffices to say that it will, as usual, be dynamically typed. The goal here is to suffer be able to play around with these building blocks somewhat faster, potentially being able to integrate them more easily with other languages.  Ultimately, I want to be able to write something like the following: 

Now, here come the juicy questions:

  • since data types in the script will only be known at runtime (i.e., once the script is being executed), how can an hypothetical interpreter going over myscript.xyz dispatch the calls to the correct overloads in the C++ functor? And,
  • how can I make a functor’s operator()s invocable in the most type-safe way possible? 

Ideally, I’d like to have something akin to Qt‘s Q_INVOKABLE tag, but I most certainly don’t want to write a Meta-Object Compiler from scratch, and I’d like to avoid devilish macros as much as possible. So, strap in for some weird templated exploration.

Boundary conditions

To better understand where this is headed, let’s define some constraints and solidify some assumptions.

I want an hypothetical interpreter to transparently deal with available functors, regardless of their implementation specifics. Furthermore, I might want to load/unload functions dynamically, e.g. via some sort of import statement. This can be solved quite classically, by having all functors inherit from the same abstract base class (like the AbstractFunctor in the above snippet), which will provide a cohesive minimal API (OOP haters intensify).

How such interpreter accesses the functors isn’t relevant here: they might be stored in a map, vector, list, etc., from which the right one is picked and called. However, there needs to be a cohesive way to pass arguments to the functor, without knowing either their types nor quantities beforehand.  
To approach this, we can take some cues from Python: when writing custom C extensions, data is passed back and forth using PyObjects, which effectively work as containers for runtime data. This is not distant from Qt’s QVariant, — also a generic runtime data container — that appears quite often in the interface between Qt/C++ and QML/JS.
To avoid rolling our own polymorphic container at first, let’s start off by using C++17’s freshly added std::variant instead.

“But”, I hear you say, “std::variant is not runtime polymorphic; all possible types need to be known at compilation time!”. That’s right. std::variant indeed defines a compile-time sum-type. However, it allows for some limited runtime introspection and provides value semantics, all which come in handy. This will temporarily restrict the data types our functors can deal with, but that’s a constraint we shall lift later. For now, let’s define a simple variant_t,

that packages a value coming from the interpreter. Let’s start off by accepting ints, doubles and strings argument types, which shall be bundled in a vector.

Our AbstractFunctor is thus bound to have a method such as void AbstractFunctor::invoke(const std::vector<variant_t>&), through which it can be, well, invoked. Notice that the functors’ return values are being ignored (our invoke method is void, after all). That’s yet another temporary simplification. We’ll get to that. Eventually.

Verifying a variant vector’s veridical validity

At this point, a first step would be asking the question: given a list of arguments (std::vector<variant_t>) and a certain function signature (operator(A, B, C, ...)), how can we verify that the argument types match said signature? 

This turned out to be an interesting point, since it lies somewhere between compile-time constraints and run-time requirements. And, while I could likely concoct some solution exploiting typeid/type_info, I figured it’d make sense to explore a path using templates. That way, I could offload as much type-checking as possible to the compilation phase. 
So, for the SumFunctor::operator(int,int) above, I’d want to write something as matches<int, int>(arglist), which verifies if arglist contains two ints. We’re thus searching for something of the form

where our typename Args... statement makes use of C++11’s template parameter packs. So, strap in for a brief jump into recursive template meta-programming territory!

Tiny trip through templated territory

Template meta-programming is basically functional programming with a cumbersome syntax. In its recursive form, it follows the principles of mathematical induction quite closely: define a base case (a proof for an n = 0) and an inductive step (proof that n = k also holds for n = k + 1), and voilà, the compiler will do the expansions for you. 

For example: in a language better suited for such task, like Haskell, we could compactly define the recursive sum of a list of integers:

Calling something like sum [1, 2, 3] would essentially expand to 1 + sum [2, 3] -> 1 + 2 + sum [3] -> 1 + 2 + 3, unsurprisingly yielding 6. To write something similar in C++, we end up being a bit more verbose:

While the syntax is unfortunately less readable, the base/inductive case pattern is very recognizable. Calling sum(1, 2, 3) spits out 6 as expected. The added constexprs aren’t mandatory, but allow for compile-time computation when possible, which is very neat.
By expanding this C++ example a bit, the sum function can be made to work with basically any combination of types for which operator+ is defined:

Sure, this introduces C++14’s decltype(auto), adds an overkill noexcept for funsies, and throws some std::move and std::forward into the mix, but hey, you get to recursively add apples and oranges together at compile time*. Tasty**!

* Provided operator+(apple, orange) is defined. YMMV.
** Could’ve been even tastier if written with a fold expression, but let’s leave that aside. 

Ok, back to our main point.

Returning to the route of recursive random ramblings

In light of this recently acquired knowledge (and after banging my head against a wall for some time), a fleshed out bool matches<...>(...) function looks like the following:

First off, there is some syntactical sugar to address here: get_if allows me to check at runtime if an underlying std::variant contains a given type. Next, operator sizeof... returns the length of a template parameter pack — we can thus know how many Types were provided as template arguments. Last, but not least, if constexpr is a C++17 godsend, that enables proper conditional branching at compile-time

In the snippet, matches acts as an entry point for our recursive checkArgumentTypeAt, which, as the name indicates, checks if the element at the N-th position of the arglist has the expected type. The base case and inductive step mechanics are the same as the defined in the previous section. Note, however, that different from our recursive sum example, the First and Second types must be explicitly peeled off from the template parameter pack. This is because the signature of checkArgumentTypeAt is the same in both specializations,  and the compiler needs to be provided with enough information for the SFINAE rules to kick in correctly. If we instead had written

since typename... Types can be empty, a call such as checkArgumentTypeAt<1, int>(arglist) would’ve been ambiguous* — as GCC happily confirms:

* Some other type deduction rules are also lurking in the shadows here. For more, here’s a mandatory Scott Meyers talk.

Deterministically dispatching to designated definitions: dope!

Once the matches function verifies the argument list, I want to be able to dispatch it to a corresponding method. At this point, I figured it’d make sense to define the concept of a Dispatcher, i.e., an object responsible for calling a particular operator() overload on our target functor.

The Dispatcher can thus be templated according to the target functor and the signature of the operator() that will be called:

To please the OOP gods (and allow runtime polymorphism), Dispatcher inherits from an AbstractDispatcher with a minimal API. The matches and checkArgumentTypeAt from before are also present, now as const members functions.

Note that Dispatcher has a member f of type Functor. That allows storing the target functor however it is best suited: reference, pointer or RAII member (which provides support for lambdas!). Right below it is an operator(), which will basically invoke the Functor‘s operator() by correctly dereferencing f. In essence, Dispatcher::operator() is just an indirection to leverage C++’s static dispatch rules, and select the correct target Functor::operator() at compile time.

The missing piece of the puzzle is the elusive dispatch method. There are two problems here: 1) how can the N-th element in the arglist be converted to the correct type and, 2) how can Dispatcher::operator() be invoked with the right number of arguments?

Problem 1 is fairly straightforward. If the matches<...>(...) function returned true, it is known that the N-th element in the list has the same type as the N-th type in the Types... template pack. To extract it, one can leverage std::tuple in a very amusing way: 

This TypeAt expression wraps Types... in a  std::tuple, and uses std::tuple_element to fetch the type of the N-th entry — all at compile time. Now, within the Dispatcher, retrieving and converting the N-th argument in arglist is just a matter of writing std::get<TypeAt<N>>(arglist.at(N)). Neat!

A solution for problem 2 builds on top on problem 1, and our complete dispatch method takes form. Unfortunately, it ends up being way less elegant then desired:

Once again, we combine sizeof... and if constexpr to only evaluate the correct conditional branch at compilation time. Unfortunately, the total amount of supported arguments needs to be added by hand, which is rather disappointing. AFAIK, there’s no way to circumvent this (and the fact that even Qt has done something similar further convinces me of that.) If any bright soul sees a better alternative here, send a pigeon

Remaining remarks

Almost. I promise.

When developing a new functor, I’d like to define which methods will be “invocable” via the dispatcher-shenanigans conjured up in the last sections. The aforementioned mythical (but not forgotten) interpreter should see an invoke method that does all the black magic internally. So, let’s circle back to the first code snippet, and define how our AbstractFunctor could look like: 

Neat little detail: different from classes/structs, template parameter pack expansion rules for functions allow packs to appear anywhere in the list (so long as the function signature provides enough cues as to what-goes-where). This can be exploited to write a little makeInvocable helper method that doesn’t require manually specifying the Functor type.
I can then add it to the constructor of our old friend, SumFunctor:

Let’s finally put everything together:

This yields the very unsurprising — but exciting — output:

If you want to play around with this, check/modify/compile the full source-code at Wandbox (or check the GitHub Gist). Have at it!


While this was a fairly long write-up, it ended up being a fairly shallow exploration in both C++ template madness and in the scripting problem outlined in the introduction. There are still quite a few open problems: how to deal with return values in the functors? How to circumvent the limitations of the introduced variant_t? Can we add support to default arguments? The list goes on. 

While a few of those I’ve got worked out already, others are widely unresolved. In any case, I’ll save all these ramblings for a post to come yet in this decade.

Edit: I absolutely cannot NOT mention ChaiScript, which is everything this post wishes to be when it grows up: a fully developed, type-safe scripting language for C++.  Absolutely fantastic code by Jason Turner & Co. If all this nonsense peaked your interest, take a look at their GitHub!

’til next time!

*Alas, additional alliterations aren’t available anymore. 

Bare-metal? Unit-tests? Why not!

Ok, so I’ll try to keep this brief. 

Some context

Unit testing is a good thing, and you should be doing it (don’t worry: that was just me talking to myself, since I miserably fail to follow my own advice). “Well, why?”, you might ask. Simply put, it’s a nice way of ensuring that your software modules/classes/thingies work in the way they should in isolation. Even though coming up with tests takes time, overall this practice tends reduces development complexity and headaches, forces you to write marginally better APIs and serves as a makeshift documentation for your code. 

Stolen from Reddit.

Good. When we’re talking about unit testing for desktop targets, framework options are plentiful. For the C/C++ world, three options pop in mind fairly quickly, but the list is gigantic. Same goes for most of the popular languages that kids use these days. 

However, in the embedded scene, the scenario is different. Computing power tends to be much more limited, and there are often interactions with external hardware that make fully isolated tests hard or unfeasible. If you’re running a *nix box, you can perhaps still profit from the one of the aforementioned frameworks (or try some tailor-made tool), but if you’re dealing with bare-metal targets, you’re pretty much out of luck. On her book, “Making Embedded Systems”, Elicia White skims over some testing practices for embedded systems, but the gist is that there’s no one-size-fits all approach. This thread pretty much echoes this feeling. 

My 2 cents

Getting to the point. As mentioned in another post, I’ve been working with a colleague on a firmware for a VFD. The firmware quickly racked quite a few software components, and keeping track of what worked and what didn’t was getting way to confusing. 

So, I drew some (*ahem*, a lot of) inspiration from Martin Dørum’s unit-test framework Snow – check it out! -, and started working on a bare-bones version of it for the STM32F10x I was using. The result is L_TEST, my attempt at an minimalist unit-test framework for bare-metal targets. The idea behind L_TEST is allowing the user to quickly define test routines and use basic assertion macros to ensure everything is working as it should. Test cases should look straightforward and clear. I really enjoy Dørum’s take on this, so L_TEST employs pretty much the same syntax as Snow does. Take a look at a simple example:

L_TEST uses printf, so you’ll need to implement the _write syscall. The output is color formatted, and should play nice with bash & friends:

Contrary to Dørum’s Snow, L_TEST is unfortunately not header-only, but the implementation is fairly light: 260 bytes for the aforementioned STM32F10x compiled with the -g flag, as shown in the map file below:

The addition of such functions also allowed the macros to be generally simpler, which in turn also reduces code footprint. Still a win, even though you have another source file to add to your list. Speaking of macros, the list still isn’t very plentiful, but here are a few:

  • L_TEST_MODULE(mname, mdesc): Define a test module, with a name and a string description. Test cases should be defined within its body
  • L_TEST_CASE(desc): Define a test case with a description. Test cases can only be defined withen a L_TEST_MODULE body.
  • L_TEST_ASSERT(val): Asserts that the value is non-zero. If the assertion fails, an error message is printed and the test case is aborted.
  • L_TEST_ASSERTEQ_INT(val, ref): Asserts that val == ref, assuming integer values. If the assertion fails, the expected and current values will be printed, and the test case will be aborted.
  • L_TEST_ASSERTEQ_FLT(val, ref): Asserts that the RPD between val and ref is less than 0.1%. If the assertion fails, the expected and current values will be printed, and the test case will be aborted.
  • L_TEST_ASSERTEQ_STR(val, ref): Executes a strcmp between the two supplied strings. If the assertion fails, the expected string and the supplied string are printed.

You can take a look at the full list here.

And thus, 

Unit-testing functionalities like this are just the tip of the iceberg. There’s a lot more that can be done using dummies, mocks, stubs and so on. L_TEST is not, by any stretch of the imagination, an attempt at a full-fledged test framework. However, it did suit my project’s simple needs, so it might suit someone else’s.  You can take a look at L_TEST on GitHub; have your go with it, and let me if you find (or better, when you find) any bugs/issues. Suggestions are welcome! 

’til next time. 


AS5030 magnetic encoder: capturing a PWM signal with an ATSAMD21

It seems that I can’t avoid periodically ostracizing this page. Welp, let’s try to make it up.

Context: the AS5030 magnetic encoder IC

In a project I’m currently working on (more about it in later posts, perhaps), I needed a halfway decent way of measuring the angular displacement of a small, manually-turnt wheel. I had originally mounted a small mechanical encoder (Bourns PEC11R) into my solution, which yielded about 96 pulses per-revolution (PPR). However, on a long run, contacting encoders are not the best pick for these applications due to mechanical wear.

Searching for a better solution I came across ams’ portfolio of magnetic encoders, which features interesting solutions for contactless position measurement. Amongst them, the AS5030 (datasheet) was the one that best met my requirements, providing an improved 256 PPR. The gist of it is simple: mount the IC, power it up and spin a magnet on top of it. The IC will generate a PWM signal with a pulse width proportional to the magnet’s angle in relation to the chip. Spoiler alert: not any magnet will work! You’ll have to get diametric magnets, such as this little guy here.

AS5030 overview: the hall frontend of the chip measures the orientation of the field lines of a diametric magnet and produces both PWM and analog signals accordingly. A SPI interface is also available.

As you may see in the diagram above, the AS5030 also provides a serial interface. The thing, however, is a weird half-duplex SPI that operates on 21-bit long frames, making it odd to use. On top of that, the AS5030 is strictly 5V compatible, which means I would have to level-shift all the signals going to my 3v3 microcontroller. That would just add unnecessary parts and complexity, so I ditched the serial interface and went with PWM (and a resistive divider to shift the voltage for the uC).

The PWM signal frequency is rated at 1.72kHz, but may vary slightly with temperature. The duty cycle encodes the position, going from a spec’d minimum of 2.26us to a spec’d maximum of 578.56us, like shown in the picture below:

AS5030 PWM signal specifictaion.

Getting the angular position accurately can be done via the ratio of the duty cycle (t_{\text{on}}) and the PWM’s period (t_{\text{on}} + t_{\text{off}}), as per the equation supplied in the datasheet:

\text{angle}[{^\circ}] = \frac{360}{256}\large[\large(257\frac{t_{on}}{t_{\text{on}}+t_{\text{off}}}\large)-1\large]

 ATSAMD21: peripheral bureaucracy

Measuring the pulse-width of a PWM signal is a textbook example of input capture. Capture operations allow you to record a signal edge together with a timestamp directly via hardware, without having to employ any CPU cycles. The vast majority of modern microcontrollers support this feature, and the ATSAMD21E18A employed in this project is no exception.

As we see in the ATSAMD21 datasheet (section 31.2), the TCC peripheral not only supports this sort of operation, but also has a dedicated mode for pulse-width capture, where “period T will be captured into CC0 and the pulse-width t_p into CC1″. Unfortunately, there’s a catch. Usually, capture operations occur entirely within a timer peripheral, which detects the signal’s edge on a dedicated pin an stores the timestamp based on its internal counter. However, to increase flexibility, capture operations in the ATSAMD21 use the External Interrupt (EIC) peripheral to generate an event in the internal Event System (EVSYS), which gets relayed via the internal event channels to the Timer/Counter (TCC) peripheral, in turn triggering the capture. Wait, what? Let me sketch that out:

Signal flow on ATSAMD21 capture operation: the PWM signal is routed from the GPIO pin to the EIC peripheral, generating an EVENT (IRQ-esque feature) on high logic. The event edges get relayed via an EVSYS channel, where the EIC acts as a generator, and the TCC timer peripheral is the event user.

This certainly allows for a lot more flexibility, not tying the capture operation to a particular pin, since all pins have EIC Interrupt/Event capability. However, configuring it is quite a mouthful (and quite poorly documented, BTW).

Talk is chip, show me the encode

I suck at puns. Ok, so let’s get a bloatw…I mean, ASF-free setup going. In my particular example, I’ll be using pin PA11, tied to EIC channel 4, but again, any pin can be used. Also, I’ll be using EVSYS’s channel 0, but you can use any of the available channels for any event. This will be done in four steps:

We’ll first configure the timer. It’ll run off the main GCLK (in this case, at 48MHz), counting up from 0 until 0xFFFFFFFF (then wrapping around):

Then, let’s setup the EIC. Sense is set to HIGH, and a tiny helper function configures the channel. Notice how the EVCTRL bit is set, which generates the EIC event:

Let’s then configure the EVSYS: EIC acts as a event generator on channel 0, and TCC1 is the event’s user (akin to a listener):

Last but not least,  let’s configure the pin. The PINMUX_PA11A_EIC_EXTINT11 value should be defined in the samd21.h header:

That’s it! Your device is now running PWM pulse-width capture on pin PA11, with no CPU cycles being used for it at all.

Now what?

Now that everything’s configured, … Well, make sure the PWM signal is getting to the designated pin. The signal’s pulse-width and period are now constantly being fetched and saved into the CC0 and CC1 registers, respectively. To wrap it up with the AS5030, we can now compute the measured angle with the following function:

This function returns angles x10, (i.e. 15.4° = 154), so that no float support needs to be added (which can eat a lot of Flash on smaller devices). It’s also enough to properly deal with the AS5030’s 1.4° resolution.

Best part? It actually works:

Rotating a magnet in front of the AS5030. Very professional test setup.

You can optionaly use the TCC_INTFLAG_MC1 interruption if you don’t want to poll the registers for changes. Yeah.

’til next time.

Getting OpenProject up and running

Ok, so I was looking for a self-hosted project management tool. Something that would fill the gap that the late (and discontinued?) DotProject left in my heart. After lots of Googlin’, I came across OpenProject. Though they offer hosting plans, you can host the tool for free on your own server. Installation looked like it had a bit too many steps for my own taste, but I went with it. Needless to say, if I’m writing this note-to-self, sh*t went south. Considering my overall stupidity regarding all-things-web, I found the installation process fairly convoluted (and apparently I’m not alone). The steps I took to get it working are sumarized below (which were compiled from quite a few different sources).

Rolling up the sleeves

A little context here – I’m installing OpenProject version 6 (latest stable release at the time of writing) on a VM running Ubuntu 14.04. Unfortunately, for other distros, YMMV. For many distros (Ubuntu included), there’s the possibility of going through with the packaged installation (e.g. via apt-get), but I ended up following the manual process. I’ll be referring to these instructions during this post. If the link breaks or changes in the future, get the PDF version (but these instructions will probably be worthless anyway):

Let’s start with the Prepare Your Environment part. I’ve totally skipped the groupadd stuff, and installed everything under my regular user – do as you please. I installed all essentials and Memcached  as listed. I already had a working install of MySQL, so I skiped the Installation of MySQL part, including the database creation (which will happen automatically later on). I then followed the Installation of Ruby (via rbenv)…

…and Installation of Node (via nodenv, but using the 6.0.0 LTS version) as listed on the manual :

Adding some swap and getting to it

I had trouble getting the rbenv install to complete. It would get killed halfway through, in which seemed an insuficcient memory issue. Anticipating a future instruction from the manual, I’ve decided to try adding some swap space to my VM. As usual, DigitalOcean has a great documentation page on how to do this. Spoilers:

Ok, now fetch the relevant version of OpenProject – I’m going with the regular openproject (not the Comunity Edition), version 6, and install it:

Now run the Configure OpenProject and Finish the Installation of OpenProject as described in the instructions. After getting them done, you’ll have an openproject/public folder that looks like this:

The lack of an index file looked super weird to me, but it’s supposed to be that way. Don’t panic.

Getting it to play nice with Apache

I started having some trouble in the Serve OpenProject with Apache section of the instructions. The compilation and installation of the passenger module worked nicely, and I could a2enmod passenger with no trouble. However, the supplied VirtualHost configuration file didn’t cut it for me.

I wanted to have OpenProject accessible from martinvb.com/op, and since my VM already serves my WordPress and Pydio, I couldn’t just move my DocumentRoot around. So I edited my  configuration in /etc/apache2/sites-enabled/ as follows:

I kept my DocumentRoot where it was, and served OpenProject through an Alias. Setting the PassengerResolveSymliksInDocumentRoot is necessary here since Passenger won’t solve symlinks automatically for versions above 2.2.0 (which happens to be the case). Also, we have to point the PassengerAppRoot to where the app’s config.ru is stored — in this case, the root of OpenProject’s cloned git repo.

Also, I’ve added rails_relative_url_root = "/op" to the config/configuration.yml created during the Configure OpenProject step, to match the folder alias I’ve created.

Reload Apache with a sudo service apache2 restart and give it a go. Works on my machine.

’til next time.

openVPS: Poor man’s motion capture

The Why and The What

Flying tiny drones indoors is cool, no questions asked. And stuff gets all the more interesting when you can accurately control the drone’s position in space — enabling all sort of crazy manouvers. However, using regular GPS for such applications is often a no-go: antennas may be too heavy for the drone, accuracy is above the often-desired cm range and getting GPS to work reliably indoors often requires breaking some laws of physics. For that reason, labs and research teams often turn to commercial marker-based motion capture systems. You place fancy cameras around a desired scene, fix some reflective makers on the objects you want to measure/track and some magic software gives you sub-centimeter (often sub-millimeter!) tracking accuracy. Yey. Too bad such systems are expensive. Like, way too expensive. Like, “I’d have to sell my organs on the black market for that”-expensive.

I thought: “surely someone wrote a piece of FOSS that allows me to grab some 3 or 4 webcams I have lying around and setup a poor-man’s version of such a system”. To my amazement, I was wrong (I mean, Google didn’t turn up anything, but I’d still love to run into something).  Considering that I still had to come up with some image-processing-related activity for my Masters degree, this seemed like a nice fit. So, let’s get on with it. [Disclaimer: heavy usage of OpenCV ahead.]

First things first: what are we looking for? Well, the sketch bellow sums it up: we have N cameras (C_0, C_1, … C_{n-1}), all looking at some region of interest from different angles. When an object of interest \{O\}  (with some red markers) shows up in the field of view of 2 \leq k \leq N cameras, we can intersect the rays emanating from each of them to discover the object’s position and orientation in space.

Sketch of system’s layout.

I’m a simple camera – I see the object, I project the object

Ok, but how? Well, let’s understand how a camera works. The image below (which I shamelessly stole) represents a pinhole-camera, a simplified mathematical camera model:

Pinhole camera model.

Any point P = (X,Y,Z) in the scene can be represented as a pixel \alpha = (u,v) via the transformation in equation (1) below. [R|t] are the rotation and translation between the camera’s frame of reference and the scene’s.

(1)   \begin{equation*} s \begin{bmatrix} u \\ v \\ 1 \end{bmatrix} = \begin{bmatrix} f_x & 0 & c_x \\ 0 & f_y & c_y \\ 0 & 0 & 1 \end{bmatrix} \bigg( R \begin{bmatrix} X \\ Y \\ Z \end{bmatrix} + t \bigg)\end{equation*}

The 3×3 matrix on the right-hand side of the equation is called the intrisic parameters matrix, and includes information about the camera’s focal distance and centering on both the and axes. As I’ve mentioned, the pinhole model is a (over)simplified one, which doesn’t accout for distortions introduced by the lens(es). Such distortions can be polynomially modelled, and the coefficients of such polynomial are then addressed as the distortion coefficients. No need to get into detail here, but both intrisic parameters and distortion coefficients can be obtained through OpenCV’s calibration routines. We’ll now assume all cameras are properly calibrated and lens distortion have been compensated, making the pinhole model valid.

Of blobs and vectors

Ok, so first thing: how to find the markers — reflective or not? Well, there’s lots of room to explore and digress here, but I’m lazy and just went with OpenCV’s built-in SimpleBlobDetector (great tutorial here!). Setup and usage are straightforward and results are nice:

Detecting markers (a.k.a. blobs).

Ok, so the detector spits out center coordinates (in pixel space, along the \{u,v\} axes) of each blob’s center point. If we want to “intersect” each view of each marker as previously sketched, then we need to compute unit vectors that emanate from each camera’s optic center \mathcal{F}_c through each blob’s center. Lets assume a point p = (x,y,z) – e.g., the marker – already in camera’s coordinates. Since the projective relationship in equation (1) loses depht information, we’ll rewrite the point in its homogenous form p_h = (x/z,  y/z, 1). We can then compute our unit vector \hat{p} of interest from the blob’s center pixel \alpha_h = (u,v,1) via,

p_h = F^{-1} \alpha_h

\hat{p} = \frac{p_h}{|p_h|}

is our aforementioned intrisic parameters matrix. Implementing this is as straightforward as it looks:

Where art thou, corresponding vector?

Following our system sketch up there, the next logical step is to intersect all these unit vectors we computed in each camera, for each marker. But the question that arises is: how do we know which vectors to match/intersect? — i.e., which vectors are pointing to the same marker on each camera? One convenient way of tackling this is exploring epipolar geometry (thanks for the inspiration, Chung J. et al). Simply put, knowing the rotation and translation between two cameras – e.g., and O’ (see image below) – allows you to map a point X as seen by O to an epiline l’ in O’. This drastically reduces the search space for the corresponding X’ — that is, the corresponding view of the same object.

There’s no need to delve into epipolar geometry — specially since I didn’t learned it in detail. OpenCV supplies satisfactory documentation and great features for this. First, I’ve used OpenCV’s stereoCalibrate within a routine that continuously checks which pairs of cameras see the calibration pattern at any given time, records that information, then later computes the fundamental and essential matrices. These which encode the [R|t] relationships between cameras, as shown in the figure below. Following that, computeCorrespondEpilines get’s you parameters a, b and c for the epiline’s equation au + bv + c = 0. So, if the n-th camera sees marker j, and a camera n+k, k \geq 1 sees a set of markers, the matching view of is m, as given by

\operatorname*{argmin}\limits_{m \in M} au_m+bv_m+c

Encoding this relations yields the results shown below (step-by-step usage of the aforementioned OpenCV functions can be found in the provided stereo_calib.cpp sample – worth checking out).

Epipolar relations between cameras being used to match points between two scenes.

Transformations between cameras

So far we have our unit vectors represented on each camera’s coordinate frame, and we are able to match the vectors “pointing” to the same marker/object on the scene. Great! Now, as a last step before computing their intersection, we need to express them all in a single coordinate frame. For simplicity, the first camera (n = 0) was chosen as the main reference for the system. So, in order to express all unit vectors in that coordinate frame, we need to know [R_n^0|t_n^0],\text{ }\forall n. In theory, we could concatenate all the successive transformations obtained with the stereo calibration procedure, but I’ve foud that it propagates lots of measurement errors and uncertainties. For now, since I’m using a small rig with only 4 cameras, the solution is to simultaneously show all cameras a known reference (namely, an AR marker), as show below:

Computing relationships between cameras via simultaneous imaging of an AR marker.

With each camera seeing the marker at the same time, we can easily compute the relationship between the 0-th camera and the n-th camera by just composing two transformations, as illustrated below:

Transformations between the AR marker, the 0-th camera and the n-th camera.

We can then easily obtain the transformations we’re interested in:

t_n^0 = -R_n^0 t_n^A + t_0^A   R_n^0 = R_0^A(R_n^A)^T

Getting to the point

Now we (finally!) have everything we need to intersect our vectors and find the position of the point(s) in the scene. I greatly recommend you take a look at this post in the Math StackExchange, since I’m piggybacking on Alberto’s detailed answer.

So, let’s assume that k \in [2, N]  cameras see a marker – thus we have unit vectors to deal with. Let’s also assume from now on, that all points and vectors are being represented in the 0-th coordinate frame. Essentially, by “intersecting” the vectors, we are trying to find a point \tilde{P} that is as close as possible to all the lines, minimizing the quadratic error relation given by

(2)   \begin{equation*}  E^2({\tilde{P}}) = \displaystyle \sum_k \bigg( ||(\tilde{P}-\mathcal{F}_{ck})||^2 - \big[ (\tilde{P} - \mathcal{F}_{ck}) \cdot \hat{p_k} \big]^2 \bigg) \end{equation*}

where \mathcal{F}_{ck} is the focal point of the k-th camera and \hat{p}_k is the unit vector emanating from the k-th camera towards our marker of interest. This, ladies and gentlemen, is also probably the most convoluted way of writing the Pythagorean Theorem. To minimize this relationship, we look for its inflection point

\frac{\delta E^2(\tilde{P})}{\delta \tilde{P}} = 0

Taking the derivative of (2) and manipulating it a bit, we get

\displaystyle \bigg[ \sum_k \big[ \hat{p_k} \cdot \hat{p_k}^T - I \big] \bigg] \cdot \tilde{P} = \displaystyle \sum_k \big[ \hat{p_k} \cdot \hat{p_k} \big] \cdot \mathcal{F}_{ck}

This can be represented as the matrix system S \tilde{P} = C. Thus our intersection point is given by \tilde{P} = S^{+} C, where S^{+} = S^T(SS^T)^{-1} is the pseudoinverse of S.

But will it blend?

How disappointing would it be, not to have a single video after all this boring maths? There you go.

If you like looking at badly-written, hastily-put-together code, there is a GitHub page for this mess. I plan on making this a nice, GUI-based thing, but who knows when. For more info, there’s also a report in portuguese, here.

’til next time.

libFilter add-ons

Following my last post on my minimalist filter library, I just got off my butt to add some high-pass filtering capabilities too. That’s specially useful when you’re trying to remove trends from datasets – this happens a lot for instance in biomedical applications, e.g. ECG, where some breathing artifacts come up as low frequency trends on your signal.

A simple example of such trend removal is included. As shown in the snippet below, taken from the example, after instancing the Filter, extracting the trend and DC components is very straigthforward:

The result looks like the following (blue senoidal wave on top is signal + trend + DC, red senoid is signal + DC, orange senoid is the filtered signal):

The filters implemented go only up to order 2 (I was lazy, sorry), and are also based on normalized Butterworth polynomials, discretized via Bilinear Transforms. Frequency responses for two hypotetical filters at 1kHz sampling frequency are shown below:

All filters are implemented in closed forms, so they can be reparameterized on the fly with minimum computational effort. Enjoy (and report any bugs). Major cleanup on the code is foreseen (and a notch-filtering example too).

’til next time.

Minimalist low-pass filter library

So, the other day I needed to compute some low-pass filters on the fly on an Arduino Mega. By “on the fly” I mean that the filters’ parameters would eventually be recomputed mid-operation, so setting some equation with static const parameters would not cut it. Using a static filter, is, however, the most common application scenario. If that’s your case – and I insist – tweak the sh*t out of your filter with a decent tool then just add it to your code. If, on the other hand, you need to update your filter dynamically (or if you’re just plain lazy to compute poles and zeros), then this is for you.

I ended up putting together a minimalist library for that, libFilter (github!). Works with the Arduino IDE and as a regular C++ library (just uncommenting one #define does the trick).

Using two filter instances on a signal from a load cell getting punched.

Using two filter instances on a signal from a load cell getting punched.

For now it only implements low-pass filters based on  normalized butterworth polynomials, but who knows what necessity might add to it next. Let’s take a look at some ADC filtering:

’til next time!

P.S.: I just can’t avoid letting this page eventually fall into temporary ostracism. Geez.

PCL library with Kinect under OSX 10.11

This last week, I dug up my trustworthy Kinect for a spin. I’ve been wanting to mess around with the PCL (Point Clouds) library for some time, so I decided to give it a shot.

Installation on OSX using Homebrew is fairly straigthforward, as shown in their documentation. However, I want to make sure that I have support for the Kinect (the Xbox 360 model).

Side note on Kinect support: To get data off of your Kinect, you can use the OpenNI library (which handles Kinect “1”).  OpenNI2 does exist, but it handles only the Kinect “2” and Occipital’s Structure Sensor. I’ll be using OpenNI here, because it’s supported directly by PCL. However, for standalone applications, I’d greatly recommend libfreenect (also available through Homebrew), which is fast, lightweight, and very easy to use.

So, we’ll be using Homebrew to get the following libraries: Boost, Flann, Qt4, GLEW, VTK and last but not least, OpenNI (1 and 2, just for argument’s sakes).  Run:

brew update
brew install boost flann qt glew vtk openni openni2

Grab a coffee, ’cause this will take some time. Now, to install PCL, you may want to check the available options (brew options pcl), then install it with at least the following settings:

brew install pcl --with-openni --with-openni2

Great, if nothing intensely wrong happened midway, you should be golden. Now, at that point I was eager to get some Kinect-demo up an running, which led me to PCL’s openni-grabber page (go ahead, visit the page). They are kind enough to supply the C++ file and a CMakeLists.txt to compile it. It compiled fine, but crashed hard every time I tried to use it, spitting the message:

Assertion failure in +[NSUndoManager _endTopLevelGroupings], /Library/Caches/com.apple.xbs/Sources/Foundation/Foundation-1256.1/Misc.subproj/NSUndoManager.m:359
2016-01-06 12:50:51.779 openni_grabber[42988:7013866] +[NSUndoManager(NSInternal) _endTopLevelGroupings] is only safe to invoke on the main thread.

Now, I don’t know if this happens under other platforms. After wasting some days struggling with this, I realized what the heck was going on. Let’s take a look at the non-working supplied code. In the header section, we see:

It uses the PCL’s Cloud Viewer, which has been long deprecated (meaning that this example is very out-of-date). Secondly, we see that the grabber’s functionality is based on a callback, defined here:

The cloud_cb_ callback function is called every time a new data packet arrives from the Kinect. This is fine, but the showCloud() command updates the display of the point cloud from within the callback, witch is what’s creating our "...is only safe to invoke on the main thread..." error. To fix that, the callback should only update the cloud’s data, while the cloud itself must be displayed with a call placed on the main thread. The fixed code I came up with looks like this:

Now, the callback cloud_cb_ only updates the cloud’s data, not touching the UI. The drawing happens in a loop inside the SimpleOpenNIViewer::run(). Since the callback and the loop are handling the same data set asynchronously, I’ve added a mutex around the critical sections to avoid any issues (the nice CPP wrapper  stolen from libfreenect’s examples, thanks!). To compile it, the following CMakeLists.txt should do the trick:

The results are not astounding: the RGB data is slightly out of place with the depth data, but I believe this is a problem internal to PCL, since my Kinect is properly calibrated. To use only the depth data in the code above, simply replace all instances of pcl::PointXYZRGBA with pcl::PointXYZ.

’til next time!

Using CMake and Qt Creator 5.5.1

Well, things have been dead around here. So, to keep things running, I’ve decided to post some less important content, mostly as notes-to-self (ya know, when you spend the weekend trying to get something to work, only to forget how you did it a month later). To avoid purging hard earned pseudo-knowledge, I’ll try to create the habit to write it down here.

So, for my first post in this (otherwise endless) series: how to use CMake with and IDE.

CMake is a nice cross-platform utility that generates Makefiles for C/C++ applications. It will automatically search for include files, libraries, generate configurable source code, etc, etc. A simple CMakeLists.txt file looks like this:

This neatly handles dependencies and configurations in a much more platform-independent way. Nice tutorials on CMake are available here.

As the kind of guy that can’t cope with VI(M) (sorry, I’m just to dumb for it), I generally gravitate towards programming-with-IDEs. I began searching for an IDE that had decent CMake integration. “Eclipse can handle it, for sure”, I thought. Nope. I was disappointed to find out that Eclipse not only lacks native CMake support, but there seems to be no plugins that handle it. CMakeBuilder was widely referenced online, but it seems that it is no longer available (using Eclipse to  query www.cmakebuilder.com/update/ yields no results).

Qt Creator ended up being the go-to solution. AFAIK, is the only well supported cross-platform GUI that handles CMake natively. Suppose we have the test_project folder, with a main.cpp and the CMakeLists.txt I’ve shown above. Use Qt Creator to open a file or project:

open_file_qt_creatorOpen the CMakeLists.txt. You’ll be prompted to run CMake (which you can skip, if you wish). After running the wizard, the left project pane will look something like this:


Neither the project name nor source files will be displayed, which is a bummer. This seems to be an issue with Qt Creator >5.0. There’s an fairly easy fix. Add the aux_source_directory command to your CMakeLists.txt:

Run Qt’s import wizard now, and the result will be as expected: project name will be featured, together with listed source files.



Great. Hit ⌘B and watch the build magic happen.

’til next time.