diff --git a/prob_laplace.m b/prob_laplace.m
index 019987eeb5f9a8a4b48b7ea0902aa835be765298..45d2ff8f6e1292103284d777e8314d0d8da3326d 100644
--- a/prob_laplace.m
+++ b/prob_laplace.m
@@ -490,43 +490,55 @@ bp = 5;
 direc = -2; % 1: jump from left to right, -1: jump from right to left.
 
 %% 
-[ls, p_i_i] = frap_distro(t, params, 2, bp, 1);
-[~, p_o_o] = frap_distro(t, params, -2, bp, 1);
-[~, p_i_o] = frap_distro(t, params, 1, bp, 1);
-[~, p_o_i] = frap_distro(t, params, -1, bp, 1);
+[ls_i_i, p_i_i] = frap_distro(t, params, 2, bp, 1);
+[ls_o_o, p_o_o] = frap_distro(t, params, -2, bp, 1);
+[ls_i_o, p_i_o] = frap_distro(t, params, 1, bp, 1);
+[ls_o_i, p_o_i] = frap_distro(t, params, -1, bp, 1);
 %%
 p_i_i(isnan(p_i_i)) = 0;
 p_o_o(isnan(p_o_o)) = 0;
 p_i_o(isnan(p_i_o)) = 0;
 p_o_i(isnan(p_o_i)) = 0;
 %%
+i = 2;
 figure; hold on;
-plot(ls, sum(p_i_i, 2));
-plot(ls, sum(p_o_o, 2));
-plot(ls, sum(p_i_o, 2));
-plot(ls, sum(p_o_i, 2));
-save('/data/biophys/hubatsch/MatlabWorkspaces/prob_laplace_190121.mat');
+plot(ls_i_i, p_i_i(:, i));
+plot(ls_o_o, p_o_o(:, i));
+plot(ls, p_i_o(:, i));
+plot(ls, p_o_i(:, i));
+%%
+figure; hold on;
+plot(ls_i_i, sum(p_i_i, 2));
+plot(ls_o_o, sum(p_o_o, 2));
+plot(ls_i_o, sum(p_i_o, 2));
+plot(ls_o_i, sum(p_o_i, 2));
+% save('/data/biophys/hubatsch/MatlabWorkspaces/prob_laplace_190121_2.mat');
 %% Save time points for comparison with python
 csvwrite('lb_in_out.csv', ls)
 csvwrite('prob_in_out.csv', p);
 %%
-KL = zeros(1, 100);
-for i = 1:68
-i_t = 1+i;
-j_t = 2+i;
 integrand = @(p, q) abs(p.*log(p./q));
 integrate = @(x, y) nansum(diff(x)'.*(y(2:end)+y(1:end-1))/2);
-n_i = @(p1, p2, p3, p4) abs(integrate(ls, p1) + integrate(ls, p2)+...
-                            integrate(ls, p3) + integrate(ls, p4));
+KL = zeros(1, 50);
+for i = 1:42
+i_t = 1+i;
+j_t = 2+i;
+n_i = @(ii, oo, io, oi) abs(integrate(ls_i_i, ii) + integrate(ls_o_o, oo)+...
+                            integrate(fliplr(ls_i_o), flipud(io)) +...
+                            integrate(ls_o_i, oi));
 N_i = n_i(p_i_i(:, i_t), p_o_o(:, i_t), p_i_o(:, i_t), p_o_i(:, i_t));
 N_j = n_i(p_i_i(:, j_t), p_o_o(:, j_t), p_i_o(:, j_t), p_o_i(:, j_t));
-
-KL(i) = integrate(ls, integrand(p_i_i(:, i_t)/N_i, p_i_i(:, j_t)/N_j))+...
-     integrate(ls, integrand(p_o_o(:, i_t)/N_i, p_o_o(:, j_t)/N_j))+...
-     integrate(ls, integrand(p_i_o(:, i_t)/N_i, p_o_i(:, j_t)/N_j))+...
-     integrate(ls, integrand(p_o_i(:, i_t)/N_i, p_i_o(:, j_t)/N_j));
+i_i = integrate(ls_i_i, integrand(p_i_i(:, i_t)/N_i,...
+                flipud(p_i_i(:, j_t))/N_j));
+o_o = integrate(ls_o_o, integrand(p_o_o(:, i_t)/N_i,...
+                flipud(p_o_o(:, j_t))/N_j));
+i_o = integrate(fliplr(ls_i_o), integrand(flipud(p_i_o(:, i_t))/...
+                N_i, p_o_i(:, j_t)/N_j));
+o_i = integrate(ls_o_i, integrand(p_o_i(:, i_t)/N_i, p_i_o(:, j_t)/N_j));
+KL(i) = i_i + o_o + i_o + o_i;
 end
-plot(t(2:101), abs(KL)/0.02);
+
+plot(t(2:51), KL/0.02);
 % ylim([0, 1])
 %%
 figure; hold on;
@@ -537,7 +549,9 @@ plot(ls, p_o_i(:, j_t)/N_j, 'b');
 %%
 plot(T_mov.t, s_dot)
 hold on; xlim([0, 2]);
-plot(t(2:101), 0.1908/0.1885*0.61*abs(KL)/0.02);
+plot(t(2:101), 0.617*abs(KL)/0.02);
+%%
+plot(t(2:51), 0.557*abs(KL)/0.04);
 %% %%%%%%%%%%%%%%%%%%% INTEGRATION FUNCTION DEFINITIONS %%%%%%%%%%%%%%%%%%%
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
@@ -558,7 +572,10 @@ end
 T_mov = Ternary_model(0, 'FRAP', params, t, 0.2);
 T_mov.solve_tern_frap();
 ls = -abs(direc)/direc*(0.001:0.04:4);
-n_T=70;
+if abs(direc) == 2 % jumps within same phase
+    ls = sort([ls, -ls]); % need pos. and neg. l, since direction not set
+end
+n_T=45;
 p = nan(length(ls), n_T);
 for j = 1:n_T
     tic