Cycles generates distorted depth
Thanks to the Blender dev forum, I figured it out.

The Cycles camera is a pinhole model that "uses the distance between a given point and the pinhole as its Z depth".
So we cannot apply the depth to points on the camera plane, but rather points that lie on the unit sphere, hence normalized.

Here is the correct formula to unproject a camera pixel into a 3D point :

    import math

    # coordinates on camera plane
    x = (0.5 - float(u) / float(img.shape[1])) * camdata.sensor_width / camdata.lens;
    y = (0.5 - float(v) / float(img.shape[0])) * sensor_height / camdata.lens;
    z = 1.0;

    norm = math.sqrt(x*x + y*y + z*z)

    # normalize = project point on unit sphere, then apply depth
    x = d * x / norm
    y = d * y / norm
    z = d * z / norm


The relevant code is in kernel_camera.h