Arbitrary camera, part 1

cuburn/affine.py

@@ -0,0 +1,34 @@
+"""
+Some simple operations on 2D affine matrices. These matrices are all stored
+in row-major order, like C, instead of Fortran-style column-major storage.
+"""
+
+import numpy as np
+
+def from_flam3(a):
+    """Convert from flam3-format [3][2] arrays to an affine matrix."""
+    return np.matrix([ [a[0][0], a[1][0], a[2][0]]
+                     , [a[0][1], a[1][1], a[2][1]]
+                     , [0, 0, 1]])
+
+def scale(x, y):
+    return np.matrix([[x,0,0], [0,y,0], [0,0,1]])
+
+def translate(x, y):
+    return np.matrix([[1,0,x], [0,1,y], [0,0,1]])
+
+def rotOrigin(rad):
+    c = np.cos(rad)
+    s = np.sin(rad)
+    return np.matrix([[c, -s, 0], [s, c, 0], [0, 0, 1]])
+
+def rotate(rad, x, y):
+    """Rotates around the given point (x, y)."""
+    return translate(x, y) * rotOrigin(rad) * translate(-x, -y)
+
+def apply(m, x, y):
+    """Apply matrix to point, returning new point as a tuple. Extends point
+    to homogeneous coordinates before applying. Mostly here as an example."""
+    r = m * np.matrix([x, y, 1]).T
+    return r[0], r[1]
+

cuburn/code/iter.py

@@ -39,7 +39,7 @@ __device__
 void apply_xf{{xfid}}(float *ix, float *iy, float *icolor,
                       const iter_info *info, mwc_st *rctx) {
     float tx, ty, ox = *ix, oy = *iy;
-    {{apply_affine('ox', 'oy', 'tx', 'ty', px, 'xf.c', 'pre')}}
+    {{apply_affine_flam3('ox', 'oy', 'tx', 'ty', px, 'xf.c', 'pre')}}
 
     ox = 0;
     oy = 0;
@@ -97,7 +97,20 @@ void iter(mwc_st *msts, iter_info *infos, float *accbuf, float *denbuf) {
 
         nsamps--;
 
-        if (x <= -0.5f || x >= 0.5f || y <= -0.5f || y >= 0.5f) {
+        // TODO: this may not optimize well, verify.
+
+        float cx, cy;
+        {{apply_affine('x', 'y', 'cx', 'cy', packer,
+                       'cp.camera_transform', 'cam')}}
+
+        float ditherwidth = {{packer.get('0.5 * cp.spatial_filter_radius')}};
+        float ditherx = mwc_next_11(&rctx) * ditherwidth;
+        float dithery = mwc_next_11(&rctx) * ditherwidth;
+
+        int ix = trunca(cx+ditherx), iy = trunca(cy+dithery);
+
+        if (ix < 0 || ix >= {{features.width}} ||
+            iy < 0 || iy >= {{features.height}} ) {
             consec_bad++;
             if (consec_bad > {{features.max_oob}}) {
                 x = mwc_next_11(&rctx);
@@ -108,11 +121,7 @@ void iter(mwc_st *msts, iter_info *infos, float *accbuf, float *denbuf) {
             continue;
         }
 
-        float ditherx = mwc_next_11(&rctx) * 0.5f;
-        float dithery = mwc_next_11(&rctx) * 0.5f;
-
-        int i = ((int)((y + 0.5f) * 1022.0f + ditherx) * 1024)
-              +  (int)((x + 0.5f) * 1022.0f + dithery) + 1025;
+        int i = iy * {{features.height}} + ix;
 
         // since info was declared const, C++ barfs unless it's loaded first
         float cp_step_frac = {{packer.get('cp_step_frac')}};
@@ -127,8 +136,8 @@ void iter(mwc_st *msts, iter_info *infos, float *accbuf, float *denbuf) {
 """)
         return tmpl.substitute(
                 features = self.features,
-                packer = self.packer.view('info'))
-
+                packer = self.packer.view('info'),
+                **globals())
 
 def render(features, cps):
     nsteps = 1000

cuburn/code/util.py

@@ -17,6 +17,18 @@ def assemble_code(*sections):
                       for kind in ['headers', 'decls', 'defs']])
 
 def apply_affine(x, y, xo, yo, packer, base_accessor, base_name):
+    return tempita.Template("""
+    {{xo}} = {{packer.get(ba + '[0,0]', bn + '_xx')}} * {{x}}
+           + {{packer.get(ba + '[0,1]', bn + '_xy')}} * {{y}}
+           + {{packer.get(ba + '[0,2]', bn + '_xo')}};
+    {{yo}} = {{packer.get(ba + '[1,0]', bn + '_yx')}} * {{x}}
+           + {{packer.get(ba + '[1,1]', bn + '_yy')}} * {{y}}
+           + {{packer.get(ba + '[1,2]', bn + '_yo')}};
+    """).substitute(x=x, y=y, xo=xo, yo=yo, packer=packer,
+                    ba=base_accessor, bn=base_name)
+
+def apply_affine_flam3(x, y, xo, yo, packer, base_accessor, base_name):
+    """Read an affine transformation in *flam3 order* and apply it."""
     return tempita.Template("""
     {{xo}} = {{packer.get(ba + '[0][0]', bn + '_xx')}} * {{x}}
            + {{packer.get(ba + '[1][0]', bn + '_xy')}} * {{y}}
@@ -42,6 +54,14 @@ uint32_t gtid() {
                 (threadIdx.z + blockDim.z *
                     (blockIdx.x + (gridDim.x * blockIdx.y))));
 }
+
+__device__
+int trunca(float f) {
+    // truncate as used in address calculations
+    int ret;
+    asm("cvt.rni.s32.f32    %0,     %1;" : "=r"(ret) : "f"(f));
+    return ret;
+}
 """
 
 class DataPackerView(object):

cuburn/render.py

@@ -9,10 +9,9 @@ from fr0stlib import pyflam3
 from fr0stlib.pyflam3._flam3 import *
 from fr0stlib.pyflam3.constants import *
 
+from cuburn import affine
 from cuburn.variations import Variations
 
-Point = lambda x, y: np.array([x, y], dtype=np.double)
-
 class Genome(pyflam3.Genome):
     @classmethod
     def from_string(cls, *args, **kwargs):
@@ -26,6 +25,25 @@ class Genome(pyflam3.Genome):
         dens /= np.sum(dens)
         self.norm_density = [np.sum(dens[:i+1]) for i in range(len(dens))]
 
+    scale = property(lambda cp: 2.0 ** cp.zoom)
+    adj_density = property(lambda cp: cp.sample_density * (cp.scale ** 2))
+    ppu = property(lambda cp: cp.pixels_per_unit * cp.scale)
+
+    @property
+    def camera_transform(cp):
+        """
+        An affine matrix which will transform IFS coordinates to image width
+        and height. Assumes that width and height are constant.
+        """
+        # TODO: when reading as a property during packing, this may be
+        # calculated 6 times instead of 1
+        return ( affine.translate(0.5 * cp.width, 0.5 * cp.height)
+               * affine.scale(cp.ppu, cp.ppu)
+               * affine.translate(-cp._center[0], -cp._center[1])
+               * affine.rotate(cp.rotate * 2 * np.pi / 360,
+                               cp.rot_center[0],
+                               cp.rot_center[1]) )
+
 class Animation(object):
     """
     Control structure for rendering a series of frames.
@@ -86,6 +104,9 @@ class Features(object):
         if any(lambda cp: cp.final_xform_enable):
             raise NotImplementedError("Final xform")
 
+        self.width = genomes[0].width
+        self.height = genomes[0].height
+
 class XFormFeatures(object):
     def __init__(self, xforms, xform_id):
         self.id = xform_id
@@ -96,27 +117,4 @@ class XFormFeatures(object):
             self.vars = (
                 self.vars.union(set([i for i, v in enumerate(x.var) if v])))
 
-class Camera(object):
-    """Viewport and exposure."""
-    def __init__(self, frame, cp):
-        # Calculate the conversion matrix between the IFS space (xform
-        # coordinates) and the sampling lattice (bucket addresses)
-        # TODO: test this code (against compute_camera?)
-        scale = 2.0 ** cp.zoom
-        self.sample_density = cp.sample_density * scale * scale
-
-        center = Point(cp._center[0], cp._center[1])
-        size = Point(cp.width, cp.height)
-
-        # pix per unit, where 'unit' is '1.0' in IFS space
-        self.ppu = Point(
-            cp.pixels_per_unit * scale / frame.pixel_aspect_ratio,
-            cp.pixels_per_unit * scale)
-        cornerLL = center - (size / (2 * self.ppu))
-        self.lower_bounds = cornerLL - gutter
-        self.upper_bounds = cornerLL + (size / self.ppu) + gutter
-        self.norm_scale = 1.0 / (self.upper_bounds - self.lower_bounds)
-        self.norm_offset = -self.norm_scale * self.lower_bounds
-        self.idx_scale = size * self.norm_scale
-        self.idx_offset = size * self.norm_offset
 

main.py

@@ -49,8 +49,6 @@ def main(args):
     noalpha = np.delete(accum, 3, axis=2)
     scipy.misc.imsave('rendered.png', noalpha)
 
-
-
     if '-g' not in args:
         return