diff --git a/scripts/run.jl b/scripts/run.jl
index 126e2b9f1d33f522c654a136d3390b10e1482a39..11668830a146273bcdc57450e379a218733fca6b 100644
--- a/scripts/run.jl
+++ b/scripts/run.jl
@@ -1,8 +1,9 @@
 include("run_experiments.jl")
 
 const datapath = joinpath(@__DIR__, "..", "data")
+const outpath = "/home/stev47/stuff/phd/data"
 
-ctx = Util.Context(datapath)
+ctx = Util.Context(datapath, outpath)
 
 ctx(experiment_convergence_rate, "fem/convergence/rate")
 ctx(experiment_optflow_middlebury_all, "fem/optflow/middlebury")
diff --git a/scripts/run_experiments.jl b/scripts/run_experiments.jl
index 0e43dace9b062a3cea2ef95a6dd967f6780f7012..170bed8166f77a0eee807aff131abbe1e06aed45 100644
--- a/scripts/run_experiments.jl
+++ b/scripts/run_experiments.jl
@@ -770,13 +770,10 @@ function denoise(img; params...)
     return st
 end
 
-
-function denoise_pd(st, img; df=nothing, algorithm, params_...)
-    params = NamedTuple(params_)
+function denoise_pd(ctx)
+    params = ctx.params
     m = 1
-    img = from_img(img) # coord flip
-    mesh = init_grid(img)
-    #mesh = init_grid(img, 5, 5)
+    mesh = ctx.params.mesh
 
     T(tdata, u) = u
     S(u, nablau) = u
@@ -795,8 +792,8 @@ function denoise_pd(st, img; df=nothing, algorithm, params_...)
     sigma = inv(tau * L^2)
     theta = 1.
 
-    #project_l2_lagrange!(st.g, img)
-    interpolate!(st.g, x -> evaluate_bilinear(img, x))
+    #project_l2_lagrange!(st.g, ctx.params.g_arr)
+    interpolate!(st.g, x -> evaluate_bilinear(ctx.params.g_arr, x))
     st.u.data .= st.g.data
 
     save_denoise(st, i) =
@@ -819,15 +816,15 @@ function denoise_pd(st, img; df=nothing, algorithm, params_...)
 
     while true
 	k += 1
-	if algorithm == :pd1
+	if ctx.params.algorithm == :pd1
 	    # no step size control
 	    step_pd2!(st; sigma, tau, theta)
-	elseif algorithm == :pd2
+	elseif ctx.params.algorithm == :pd2
 	    theta = 1 / sqrt(1 + 2 * mu * tau)
 	    tau *= theta
 	    sigma /= theta
 	    step_pd2!(st; sigma, tau, theta)
-	elseif algorithm == :newton
+	elseif ctx.params.algorithm == :newton
 	    step!(st)
 	end
 	#pvd[k] = save_denoise(st, k)
@@ -836,7 +833,7 @@ function denoise_pd(st, img; df=nothing, algorithm, params_...)
 	norm_step_ = norm_step(st) * domain_factor
         #estimate_pd!(st)
 
-        log!(df; k,
+        log!(ctx.params.df; k,
             norm_step = norm_step_,
             #est = estimate_error(st),
             primal_energy = primal_energy(st))
@@ -854,24 +851,24 @@ end
 function experiment_convergence_rate(ctx)
     img = loadimg(joinpath(ctx.indir, "input.png"))
     #img = [0. 0.; 1. 0.]
+    g_arr = from_img(img)
+    mesh = init_grid(g_arr;)
 
     algparams = (
         alpha1=0., alpha2=30., lambda=1., beta=0.,
         gamma1=1e-2, gamma2=1e-3,
         tol = 1e-10,
         max_iters = 10000,
-        )
+    )
+    algctx = Util.Context(ctx; g_arr, mesh, algparams...)
 
     df1 = DataFrame()
     df2 = DataFrame()
     df3 = DataFrame()
 
-    st1 = denoise_pd(ctx, img;
-        algorithm=:pd1, df = df1, algparams...);
-    st2 = denoise_pd(ctx, img;
-        algorithm=:pd2, df = df2, algparams...);
-    st3 = denoise_pd(ctx, img;
-        algorithm=:newton, df = df3, algparams...);
+    st1 = denoise_pd(Util.Context(algctx; algorithm = :pd1, df = df1));
+    st2 = denoise_pd(Util.Context(algctx; algorithm = :pd2, df = df2));
+    st3 = denoise_pd(Util.Context(algctx; algorithm = :newton, df = df3));
 
     energy_min = min(minimum(df1.primal_energy), minimum(df2.primal_energy),
 		     minimum(df3.primal_energy))
@@ -1216,7 +1213,6 @@ function experiment_optflow_middlebury(ctx)
     maxflow = OpticalFlowUtils.maxflow(gtflow)
 
     ctx = Util.Context(ctx; imgf0, imgf1, maxflow)
-    mkpath(ctx.outdir)
     saveimg(joinpath(ctx.outdir, "ground_truth.png"), colorflow(gtflow; maxflow))
     return optflow(ctx)
 end
diff --git a/scripts/util.jl b/scripts/util.jl
index 2cc4a5b53c78a836320505f6aaa67f58d24194f5..790a51cf80fef51919c0b90d28ce2c204872f0ef 100644
--- a/scripts/util.jl
+++ b/scripts/util.jl
@@ -38,6 +38,7 @@ parameters `params`.
 function (ctx::Context)(f::Function, subdir::String; params...)
     # maybe debug output with params?
     subctx = Context(ctx, subdir; params...)
+    mkpath(subctx.outdir)
     println("starting experiment `$(subctx.path)` ...")
     f(subctx)
 end