diff --git a/scripts/run_experiments.jl b/scripts/run_experiments.jl index a8d7e39076aefada81db65aef8d1beb43bc1bda0..a59d86161a28bdd440073bb096857245e23560a2 100644 --- a/scripts/run_experiments.jl +++ b/scripts/run_experiments.jl @@ -917,12 +917,9 @@ function optflow(ctx) ctx.params.alpha1, ctx.params.alpha2, ctx.params.lambda, ctx.params.beta, ctx.params.gamma1, ctx.params.gamma2) - u_acc = FeFunction(st.u.space, name = "u_acc") - u_acc.data .= 0 - function warp!() # warp image into imgfw / fw - imgfw = warp_backwards(imgf1, sample(u_acc)) + imgfw = warp_backwards(imgf1, sample(st.u)) project_l2_pixel!(fw, imgfw) # recompute optflow operator T based on u0 and fw @@ -947,12 +944,10 @@ function optflow(ctx) project_l2_pixel!(f1, imgf1) end - interpolate_image_data!() - warp!() # just to fill st.g save_step(i) = output(st, joinpath(ctx.outdir, "output_$(lpad(i, 5, '0')).vtu"), - st.g, st.u, st.p1, st.p2, st.est, f0, f1, fw, u_acc) + st.g, st.u, st.p1, st.p2, st.est, f0, f1, fw) i = 0 pvd = paraview_collection(joinpath(ctx.outdir, "output.pvd")) do pvd @@ -962,6 +957,10 @@ function optflow(ctx) #norm_g_old = norm_l2(st.g) for j in 1:5 + println("interpolate image data ...") + interpolate_image_data!() + println("warp ...") + warp!() # in the first step just to fill st.g i += 1 # interior newton loop @@ -984,69 +983,22 @@ function optflow(ctx) #display(plot(colorflow(to_img(sample(st.u)); ctx.params.maxflow))) end - #estimate!(st) - #pvd[i] = save_step(i) - - u_acc.data .= st.u.data - display(plot(colorflow(to_img(sample(u_acc)); ctx.params.maxflow))) - - # poor man's hash - #println(integrate(st.mesh, (x; tdata) -> tdata; st.tdata)) - - #println("refine ...") - ##marked_cells = mark(st; theta = 0.5) - #marked_cells = Set(axes(mesh.cells, 2)) - #mesh, fs = refine(mesh, marked_cells; - # st.est, st.g, st.u, st.p1, st.p2, st.du, st.dp1, st.dp2, - # f0, f1, fw, u_acc) - #st = L1L2TVState(mesh, st.d, st.m, T, nabla(fs.fw), S, - # st.alpha1, st.alpha2, st.beta, st.lambda, st.gamma1, st.gamma2, - # fs.est, fs.g, fs.u, fs.p1, fs.p2, fs.du, fs.dp1, fs.dp2) - #f0, f1, fw, u_acc = (fs.f0, fs.f1, fs.fw, fs.u_acc) - - #i += 1 - - #println("interpolate image data ...") - #interpolate_image_data!() - - + display(plot(colorflow(to_img(sample(st.u)); ctx.params.maxflow))) - println("warp ...") - ##warp!() # yay estimate!(st) pvd[i] = save_step(i) - #u_acc.data .= st.u.data - warp!() # yay + println("refine ...") + marked_cells = mark(st; theta = 0.5) + #marked_cells = Set(axes(mesh.cells, 2)) + mesh, fs = refine(mesh, marked_cells; + st.est, st.g, st.u, st.p1, st.p2, st.du, st.dp1, st.dp2, + st.tdata, f0, f1, fw) + st = L1L2TVState(mesh, st.d, st.m, st.T, fs.tdata, st.S, + st.alpha1, st.alpha2, st.beta, st.lambda, st.gamma1, st.gamma2, + fs.est, fs.g, fs.u, fs.p1, fs.p2, fs.du, fs.dp1, fs.dp2) + f0, f1, fw = (fs.f0, fs.f1, fs.fw) end - - - #warp!() - #norm_g = norm_l2(st.g) - #i += 1 - #pvd[i] = save_step(i) - #i >= 50 && break - #continue - - #marked_cells = mark(st; theta = 0.5) - #marked_cells = Set(axes(mesh.cells, 2)) - - #println("refining ...") - - #mesh, fs = refine(mesh, marked_cells; - # st.est, st.g, st.u, st.p1, st.p2, st.du, st.dp1, st.dp2, - # f0, f1, fw) - #st = L1L2TVState(mesh, st.d, st.m, T, nabla(fs.fw), S, - # st.alpha1, st.alpha2, st.beta, st.lambda, st.gamma1, st.gamma2, - # fs.est, fs.g, fs.u, fs.p1, fs.p2, fs.du, fs.dp1, fs.dp2) - #f0, f1, fw = (fs.f0, fs.f1, fs.fw) - # i += 1 - # pvd[i] = save_step(i) - - #println("interpolate image data ...") - #interpolate_image_data!() - # i += 1 - # pvd[i] = save_step(i) end #display(plot(colorflow(to_img(sample(st.u)); ctx.params.maxflow)))