George MacKerron: code blog

GIS, software development, and other snippets

Lightweight inline formatting for NSAttributedString

For Cocoa developers: I’ve just put a small category on NSMutableAttributedString on Github.

It applies *bold*, /italic/, _underline_, -strikethrough-, ^superscript^ and ~subscript~ styles, and handles */nested/* and *overlapping /styles* properly/.

More details at https://github.com/jawj/NSMutableAttributedString-InlineStyles

Written by George

March 11th, 2014 at 4:13 pm

Posted in iPhone,Mac

Using vDSP_deq22 as a bandpass filter

For reasons that may become clearer in future, I needed to use a bandpass filter in an iOS app. The DSP part of Apple’s Accelerate framework makes this lightning fast both for the programmer to implement and for the machine to execute … if the programmer knows how use vDSP_deq22, which is documented, at best, concisely.

The following functions produce the five-coefficient filter definition you need to pass to vDSP_deq22. Since this took me a little while to put together, I thought I’d share.

void makeBandpassFilterWithFcAndQ(float* filter, double Fs, double Fc, double Q) {
 
  // Fs = sampling rate, Fc = centre freq
  // with thanks to http://www.earlevel.com/main/2013/10/13/biquad-calculator-v2/
  // and https://github.com/bartolsthoorn/NVDSP
 
  double K = tan(M_PI * Fc / Fs);
  double norm = 1.0 / (1.0 + K / Q + K * K);
  filter[0] = (float)(K / Q * norm);
  filter[1] = 0.0f;
  filter[2] = -filter[0];
  filter[3] = (float)(2.0 * (K * K - 1.0) * norm);
  filter[4] = (float)((1.0 - K / Q + K * K) * norm);
}
 
void makeBandpassFilterWithFreqRange(float* filter, double Fs, double Fbtm, double Ftop) {
 
  // with thanks to 
  // http://stackoverflow.com/questions/15627013/how-do-i-configure-a-bandpass-filter
  // -- this sets Q such that there's -3dB gain (= 50% power loss) at Fbtm and Ftop
 
  double Fc = sqrt(Fbtm * Ftop);
  double Q = Fc / (Ftop - Fbtm);
  makeBandpassFilterWithFcAndQ(filter, Fs, Fc, Q);
}

And you use it like so:

float filter[5];
makeBandpassFilterWithFreqRange(filter, sampleRate, filterLoRate, filterHiRate);
vDSP_deq22(rawFloats, rawStride, filter, filteredFloats, filteredStride, numRawFloats - 2); 
  // rawFloats and filteredFloats are pointers, of course

Written by George

February 4th, 2014 at 2:37 pm

Posted in iPhone,Mac

Turning a Xen domU into a VMware VM (or: How to bring your Linode home)

I recently posted a HOWTO based on my experience moving a Xen domU from Linode to my own Xen Dom0 setup at Hetzner.

Since this machine is only a development server, I more recently decided to turn the same machine into a VMware VM, running locally (in VMware Fusion 4 on my MacBook Pro). Here, I document the steps necessary for that transformation.

More »

Written by George

July 29th, 2013 at 10:26 am

Posted in System admin

Map algebra with callbacks in PostGIS 2.1

I’ve been playing with the much-enhanced raster support in PostGIS 2.1 in the last few days. Amongst other things I’ve been producing maps with ST_MapAlgebra, which has changed a little since 2.0.

The example callback function in the documentation is somewhat unhelpfully bare-bones: it always returns 0, touching none of its arguments. So here’s an equally useless — but rather more informative — example of a callback function for the simple one-raster, one-band, one-pixel-at-a-time case:

create or replace function 
callback_fn(pixel float[], pos integer[], variadic userargs text[]) 
  returns float 
  language plpgsql 
  immutable  -- careful: this function is immutable, yours may not be
as $$
  declare
    pixval float;
    inputx integer;
    inputy integer;
  begin
    pixval := pixel[1][1][1];  -- pixel indices: [raster #][xdistance][ydistance]
    inputx := pos[1][1];       -- pos indices:   [raster #][x = 1, y = 2]
    inputy := pos[1][2];       --                (raster #0 is the output raster)

    return pixval + inputx + inputy;
  end;
$$;

-- example call:

select st_mapalgebra(rast, 1, 'callback_fn(float[], integer[], text[])'::regprocedure) 
from raster_table;

And here’s something that caught me out: unless you’re passing some userargs, don’t declare your callback function to be strict. If you do, your callback function will never be called because its userargs argument is always null.

Written by George

July 21st, 2013 at 2:42 pm

Posted in GIS,PostGIS,SQL

Sunrise and sunset times with PostGIS and PL/R

Since light can affect happiness, two important pieces of environmental data I add to the Mappiness data set during analysis are: (1) whether a response was made in daylight; and (2) day length when and where the response was made.

To derive these, I need the date, time and location of the response, and sunrise and sunset times for that date and location. R’s StreamMetabolism library provides sunrise/sunset calculations based on NOAA routines. And since my data is in PostGIS, it’s handy to use PL/R to access these.

I set this up as follows on my Ubuntu 12.04 GIS box.

In bash:

sudo aptitude install postgresql-9.1-plr
sudo R

In R:

install.packages('StreamMetabolism', dependencies = TRUE)

In Postgres (9.1):

create extension plr;
create table plr_modules (modseq int4, modsrc text);
insert into plr_modules values (0, 'library("StreamMetabolism")');
 
create or replace function 
  _r_daylight_period(lat double precision, lon double precision, datestring text, 
                     timezone text)
returns setof integer as $$
  return(as.integer(sunrise.set(lat, lon, datestring, timezone)))
$$ language plr immutable strict;
 
create or replace function 
  sunrise(location geometry, moment timestamp without time zone, timezone text)
returns timestamp without time zone as $$
  select to_timestamp(min(s)) at time zone $3 from _r_daylight_period(
    st_y(st_transform($1, 4326)), -- 4326 = WGS84
    st_x(st_transform($1, 4326)),
    to_char($2, 'YYYY/MM/DD'),
    $3
  ) s
$$ language sql immutable strict;
 
create or replace function 
  sunset(location geometry, moment timestamp without time zone, timezone text)
returns timestamp without time zone as $$
  select to_timestamp(max(s)) at time zone $3 from _r_daylight_period(
    st_y(st_transform($1, 4326)), -- 4326 = WGS84
    st_x(st_transform($1, 4326)),
    to_char($2, 'YYYY/MM/DD'),
    $3
  ) s
$$ language sql immutable strict;
 
create or replace function 
  is_daylight(location geometry, moment timestamp without time zone, timezone text) 
returns boolean as $$
  select (sunrise($1, $2, $3), sunset($1, $2, $3)) overlaps ($2, cast('0' as interval));
$$ language sql immutable strict;
 
create or replace function 
  hours_in_the_day(location geometry, moment timestamp without time zone, timezone text) 
returns double precision as $$
  select extract(epoch from sunset($1, $2, $3) - sunrise($1, $2, $3)) / 3600.0;
$$ language sql immutable strict;

These new functions can be used like so:

select 
  is_daylight(geom, moment, 'Europe/London'), 
  hours_in_the_day(geom, moment, 'Europe/London') 
from mytable;

(Note: what I actually do with the Mappiness data is to use a single call to _r_daylight_period, and calculate the other quantities as a second step. This speeds things up a lot, because I have millions of rows to process and Postgres doesn’t appear to do as much caching of immutable function results as one would like here).

Written by George

October 15th, 2012 at 5:52 pm

Posted in GIS,PL/R,PostGIS,SQL

Fixing Bluetooth sleep issues with MacBook and Magic Trackpad

We have an old MacBook (the original white Intel model from 2006) running EyeTV as our telly.

Apple’s Magic Trackpad makes a handy remote control for this setup. Unfortunately, Bluetooth on the old MacBook is highly temperamental. Built-in Bluetooth regularly fails: out of the blue, the Mac decides that it has no Bluetooth module after all, puts a wavy line through the menu bar icon, and ignores the trackpad. The only fix for this is to shut the computer down, and leave it off for several minutes. This is annoying.

I had an old D-Link Bluetooth dongle (DBT-120), so I tried using this instead. This is better, in that the failure mode is less annoying. Using the dongle, the trackpad’s Bluetooth connection only fails after a prolonged sleep. It looks as if it’s still connected, with a dotted line across the Bluetooth icon, but it’s unresponsive. This can be fixed by simply unplugging and replugging the dongle. But that’s still a pain, and is rather a compromise of the ‘remote’ in remote control.

It turns out that the post-sleep unresponsiveness may also be fixed by restarting the Bluetooth daemon, by typing sudo killall blued in the Terminal.

This is good news, because we can automate this action using sleepwatcher.

If typing sudo killall blued in Terminal solves your Bluetooth issues after sleep, then you may want to use sleepwatcher too. (Note that after a sudo command, you may be asked for a password. Nothing will show up as you type, but that’s OK: just type your password and press Return).

More »

Written by George

October 11th, 2012 at 1:55 pm

Posted in Mac

Using the OSTN02 transformation in PostGIS

When converting coordinates between WGS84 (GPS) and OSGB36 (UK National Grid), using OSTN02 can gain us a few metres in accuracy over the basic parametric transformation provided by PostGIS’s st_transform via PROJ.4.

Happily, Ordnance Survey now distribute an NTv2 version of the OSTN02 transformation, courtesy of the Defence Geographic Centre, which can be used by PROJ.4 and, therefore, PostGIS.

More »

Written by George

July 3rd, 2012 at 12:59 pm

Posted in GIS,PostGIS

Objective-C SHA1 categories for NSData and NSString

I recently needed to calculate a SHA1 hash in an iOS app.

In iOS4+ it’s possible to use CommonCrypto, but Mappiness has always supported iOS3. I therefore added NSData and NSString categories to a public domain C implementation instead. This remains public domain: do with it what you will.

It relies on a hex-encoding category on NSData, which you can also consider public domain.

More »

Written by George

June 19th, 2012 at 10:12 am

Posted in iPhone

How to compile PostGIS 2.1 on Ubuntu Server 12.04+

PostGIS 2 has some exciting new goodies — including raster support — that I’m keen to use in the analysis of Mappiness data.

But the PostGIS package provided by Ubuntu 12.04 is still only at version 1.5, and the GEOS and GDAL packages are also too old to support the new version.

So — this is how I compiled PostGIS 2.0.1 2.0.3 2.1.0rc2 2.1.0 and its dependencies on my GIS server.

More »

Written by George

June 1st, 2012 at 6:30 pm

Spindlytext: write in the sky with the Google Earth API

I’ve just released Spindlytext on Github. It’s the library that powers the live data display in Pigeon Sim, by creating KML linestrings in the shape of letters.

More »

Written by George

May 2nd, 2012 at 5:30 pm