diff --git a/bin/barplot_coarse.py b/bin/barplot_coarse.py
index 67de3aefb4e9bb0d4ef42a2b61a51fd28ddc873d..47906024390a599c882437cf439a50685c45d430 100644
--- a/bin/barplot_coarse.py
+++ b/bin/barplot_coarse.py
@@ -19,19 +19,15 @@ opframe = pandas.read_csv("./operations.csv",
 
 av = opframe[opframe.what == "alpha_volume_kernel"]
 ask = opframe[opframe.what == "alpha_skeleton_kernel"]
-ab = opframe[opframe.what == "alpha_boundary_kernel"]
 deg = [i - 0.3 for i in sorted(av['degree'])]
 
 _, y1 = list(itertools.izip(*sorted(itertools.izip(av['degree'], av['ops']))))
 _, y2 = list(itertools.izip(*sorted(itertools.izip(ask['degree'], ask['ops']))))
 y2 = [a+b for a, b in zip(y1, y2)]
-_, y3 = list(itertools.izip(*sorted(itertools.izip(ab['degree'], ab['ops']))))
-y3 = [a+b for a, b in zip(y2, y3)]
-y4 = [1.0] * len(deg)
+y3 = [1.0] * len(deg)
 
-r4 = ax.bar(deg, y4, width, color="grey") 
-r3 = ax.bar(deg, y3, width, color="green") 
-r2 = ax.bar(deg, y2, width, color="blue") 
+r3 = ax.bar(deg, y3, width, color="grey")
+r2 = ax.bar(deg, y2, width, color="blue")
 r1 = ax.bar(deg, y1, width, color="red")
 
 ax.set_ylabel("Percentage")
@@ -46,19 +42,15 @@ timeframe = pandas.read_csv("./timeratios.csv",
 
 av = timeframe[timeframe.what == "alpha_volume_kernel"]
 ask = timeframe[timeframe.what == "alpha_skeleton_kernel"]
-ab = timeframe[timeframe.what == "alpha_boundary_kernel"]
 deg = [i + 0.3 for i in sorted(av['degree'])]
 
 _, y1 = list(itertools.izip(*sorted(itertools.izip(av['degree'], av['ops']))))
 _, y2 = list(itertools.izip(*sorted(itertools.izip(ask['degree'], ask['ops']))))
 y2 = [a+b for a, b in zip(y1, y2)]
-_, y3 = list(itertools.izip(*sorted(itertools.izip(ab['degree'], ab['ops']))))
-y3 = [a+b for a, b in zip(y2, y3)]
-y4 = [1.0] * len(deg)
+y3 = [1.0] * len(deg)
 
-r4 = ax.bar(deg, y4, width, color="grey", label="PDELab overhead") 
-r3 = ax.bar(deg, y3, width, color="green", label="Boundary integrals") 
-r2 = ax.bar(deg, y2, width, color="blue", label="Skeleton integrals") 
+r3 = ax.bar(deg, y3, width, color="grey", label="PDELab overhead")
+r2 = ax.bar(deg, y2, width, color="blue", label="Skeleton integrals")
 r1 = ax.bar(deg, y1, width, color="red", label="Volume integrals")
 
 ax.legend(loc=3,
diff --git a/bin/barplot_fine.py b/bin/barplot_fine.py
index 4e050c80ef775780c559ed82bfa6f4a2acba1564..47f4514da839afc43d7c588ce4a7158e5d08e5b0 100644
--- a/bin/barplot_fine.py
+++ b/bin/barplot_fine.py
@@ -28,7 +28,8 @@ av3 = opframe[opframe.what == "alpha_volume_kernel_stage3"]
 ask1 = opframe[opframe.what == "alpha_skeleton_kernel_stage1"]
 ask2 = opframe[opframe.what == "alpha_skeleton_kernel_quadratureloop"]
 ask3 = opframe[opframe.what == "alpha_skeleton_kernel_stage3"]
-ab = opframe[opframe.what == "alpha_boundary_kernel"]
+setup1 = opframe[opframe.what == "alpha_volume_kernel_setup"]
+setup2 = opframe[opframe.what == "alpha_skeleton_kernel_setup"]
 deg = [i - width/1.8 for i in sorted(av1['degree'])]
 
 
@@ -38,23 +39,26 @@ def update(frame, result=None):
         return y
     else:
         return [a+b for a, b in zip(y, result)]
-        
+
+
 y1 = update(av1)
 y2 = update(av2, y1)
 y3 = update(av3, y2)
-y4 = update(ask1, y3)
-y5 = update(ask2, y4)
-y6 = update(ask3, y5)
-y7 = update(ab, y6)
-y8 = [1.0] * len(deg)
-
-r8 = ax.bar(deg, y8, width, color="grey") 
-r7 = ax.bar(deg, y7, width, color="green") 
-r6 = ax.bar(deg, y6, width, color="xkcd:light blue") 
-r5 = ax.bar(deg, y5, width, color="xkcd:blue")
-r4 = ax.bar(deg, y4, width, color="xkcd:dark blue") 
-r3 = ax.bar(deg, y3, width, color="yellow") 
-r2 = ax.bar(deg, y2, width, color="orange") 
+y4 = update(setup1, y3)
+y5 = update(ask1, y4)
+y6 = update(ask2, y5)
+y7 = update(ask3, y6)
+y8 = update(setup2, y7)
+y9 = [1.0] * len(deg)
+
+r9 = ax.bar(deg, y9, width, color="grey")
+r8 = ax.bar(deg, y8, width, color="xkcd:light green")
+r7 = ax.bar(deg, y7, width, color="xkcd:light blue")
+r6 = ax.bar(deg, y6, width, color="xkcd:blue")
+r5 = ax.bar(deg, y5, width, color="xkcd:dark blue")
+r4 = ax.bar(deg, y4, width, color="green")
+r3 = ax.bar(deg, y3, width, color="yellow")
+r2 = ax.bar(deg, y2, width, color="orange")
 r1 = ax.bar(deg, y1, width, color="red")
 
 ax.set_ylabel("Percentage")
@@ -73,25 +77,28 @@ av3 = timeframe[timeframe.what == "alpha_volume_kernel_stage3"]
 ask1 = timeframe[timeframe.what == "alpha_skeleton_kernel_stage1"]
 ask2 = timeframe[timeframe.what == "alpha_skeleton_kernel_quadratureloop"]
 ask3 = timeframe[timeframe.what == "alpha_skeleton_kernel_stage3"]
-ab = timeframe[timeframe.what == "alpha_boundary_kernel"]
+setup1 = timeframe[timeframe.what == "alpha_volume_kernel_setup"]
+setup2 = timeframe[timeframe.what == "alpha_skeleton_kernel_setup"]
 deg = [i + width/1.8 for i in sorted(av1['degree'])]
 
 y1 = update(av1)
 y2 = update(av2, y1)
 y3 = update(av3, y2)
-y4 = update(ask1, y3)
-y5 = update(ask2, y4)
-y6 = update(ask3, y5)
-y7 = update(ab, y6)
-y8 = [1.0] * len(deg)
-
-r8 = ax.bar(deg, y8, width, label="PDELab overhead", color="grey") 
-r7 = ax.bar(deg, y7, width, label="Boundary", color="green") 
-r6 = ax.bar(deg, y6, width, label="Skeleton, backward SF", color="xkcd:light blue") 
-r5 = ax.bar(deg, y5, width, label="Skeleton, Quadrature", color="xkcd:blue")
-r4 = ax.bar(deg, y4, width, label="Skeleton, forward SF", color="xkcd:dark blue") 
-r3 = ax.bar(deg, y3, width, label="Volume, backward SF", color="yellow") 
-r2 = ax.bar(deg, y2, width, label="Volume, Quadrature", color="orange") 
+y4 = update(setup1, y3)
+y5 = update(ask1, y4)
+y6 = update(ask2, y5)
+y7 = update(ask3, y6)
+y8 = update(setup2, y7)
+y9 = [1.0] * len(deg)
+
+r9 = ax.bar(deg, y9, width, label="PDELab overhead", color="grey")
+r8 = ax.bar(deg, y8, width, label="Geometry evaluations Skeleton", color="xkcd:light green")
+r7 = ax.bar(deg, y7, width, label="Skeleton, backward SF", color="xkcd:light blue")
+r6 = ax.bar(deg, y6, width, label="Skeleton, Quadrature", color="xkcd:blue")
+r5 = ax.bar(deg, y5, width, label="Skeleton, forward SF", color="xkcd:dark blue")
+r4 = ax.bar(deg, y4, width, label="Geometry evaluations Volume", color="green")
+r3 = ax.bar(deg, y3, width, label="Volume, backward SF", color="yellow")
+r2 = ax.bar(deg, y2, width, label="Volume, Quadrature", color="orange")
 r1 = ax.bar(deg, y1, width, label="Volume, forward SF", color="red")
 
 lgd = ax.legend(loc=3,
diff --git a/bin/process_measurements.py b/bin/process_measurements.py
index 0d45c7293befca1c166617746c13069c3e9c57cc..53fb0873dee3a3cf3719c2faf5f98bc2a2925a18 100755
--- a/bin/process_measurements.py
+++ b/bin/process_measurements.py
@@ -22,10 +22,27 @@ def get_reference_kernel(kernel):
     return kernel
 
 
+def get_reference_integral_kernel(kernel):
+    if kernel.startswith("alpha") and len(kernel.split("_")) == 4:
+        return "_".join(kernel.split("_")[:-1])
+    if kernel.startswith("jacobian_apply") and len(kernel.split("_")) == 5:
+        return "_".join(kernel.split("_")[:-1])
+    return kernel
+
+
+def get_siblings_kernel(kernel):
+    s = kernel.split("_")
+    if (kernel.startswith("alpha") and len(s) == 4) or (kernel.startswith("jacobian_apply") and len(s) == 5):
+        for suffix in ["setup", "stage1", "quadratureloop", "stage3"]:
+            yield "_".join(s[:-1] + [suffix])
+    else:
+        yield kernel
+
 def calculate_operations_percentage():
-    frame = pandas.read_csv('timings.csv', header=None, names=('rank', 'ident', 'kernel', 'what', 'value'), delimiter=' ')
+    frame = pandas.read_csv('timings.csv', header=None, names=('rank', 'level', 'ident', 'kernel', 'what', 'value'), delimiter=' ')
     ops = frame[frame.what != "time"]
     ops = ops.groupby(('rank', 'ident', 'kernel'))['value'].max().to_frame().reset_index().groupby(('ident', 'kernel'))['value'].max()
+
     with open('operations.csv', 'w') as out:
         for key in ops.keys():
             ident, kernel = key
@@ -34,18 +51,26 @@ def calculate_operations_percentage():
 
 
 def calculate_times_percentage():
-    frame = pandas.read_csv('timings.csv', header=None, names=('rank', 'ident', 'kernel', 'what', 'value'), delimiter=' ')
-    time = frame[frame.what == "time"]
-    time = time.groupby(('rank', 'ident', 'kernel'))['value'].min().to_frame().reset_index().groupby(('ident', 'kernel'))['value'].max()
+    frame = pandas.read_csv('timings.csv', header=None, names=('rank', 'level', 'ident', 'kernel', 'what', 'value'), delimiter=' ')
+
+    time4 = frame[(frame.what == "time") & (frame.level == 4)]
+    time4 = time4.groupby(('rank', 'ident', 'kernel'))['value'].min().to_frame().reset_index().groupby(('ident', 'kernel'))['value'].max()
+    time3 = frame[(frame.what == "time") & (frame.level == 3)]
+    time3 = time3.groupby(('rank', 'ident', 'kernel'))['value'].min().to_frame().reset_index().groupby(('ident', 'kernel'))['value'].max()
+
     with open('timeratios.csv', 'w') as out:
-        for key in time.keys():
+        for key in time4.keys():
             ident, kernel = key
             degree = re.match(".*deg([0-9]*).*", ident).group(1)
-            out.write(" ".join([ident, degree, kernel, str(time[ident][kernel] / time[ident][get_reference_kernel(kernel)]) + "\n"]))
+            stage = time4[ident][kernel]
+            sum4 = sum(time4[ident][k] for k in get_siblings_kernel(kernel))
+            int3 = time3[ident][get_reference_integral_kernel(kernel)]
+            ref =  time3[ident][get_reference_kernel(kernel)]
+            out.write(" ".join([ident, degree, kernel, str(float(stage*int3) / float(sum4*ref)) + "\n"]))
 
 
 def calculate_floprate():
-    frame = pandas.read_csv('timings.csv', header=None, names=('rank', 'ident', 'kernel', 'what', 'value'), delimiter=' ')
+    frame = pandas.read_csv('timings.csv', header=None, names=('rank', 'level', 'ident', 'kernel', 'what', 'value'), delimiter=' ')
     time = frame[frame.what == "time"]
     ops = frame[frame.what != "time"]
 
@@ -56,11 +81,12 @@ def calculate_floprate():
         for key in time.keys():
             ident, kernel = key
             degree = re.match(".*deg([0-9]*).*", ident).group(1)
-            out.write(" ".join([ident, degree, kernel, str((ops[ident][kernel] / time[ident][kernel]) / 1e9)]) + "\n")
+            operations = ops[ident][kernel]
+            out.write(" ".join([ident, degree, kernel, str((operations / time[ident][kernel]) / 1e9)]) + "\n")
 
 
 def calculate_doftimes():
-    frame = pandas.read_csv('timings.csv', header=None, names=('rank', 'ident', 'kernel', 'what', 'value'), delimiter=' ')
+    frame = pandas.read_csv('timings.csv', header=None, names=('rank', 'level', 'ident', 'kernel', 'what', 'value'), delimiter=' ')
     dofs = frame[frame.what == "dofs"]
     time = frame[frame.what == "time"]