 }
 
 //Gabor noise based on Lagae, Lefebvre, Drettakis, Dutre 2011
-  float eval_cell(in vec2 cpos, in ivec2 gpos, in ivec2 dnbr, gnoise_params params){
+  float eval_cell(const vec2 cpos, const ivec2 gpos, const ivec2 dnbr, const gnoise_params params, const gnoise_im_params im_params){
 	float u = seed(bound_grid(gpos+dnbr)); //deterministic seed for nbr cell
-	int impulses = poisson(u, params.density); //number of impulses in nbr cell
+	int impulses = int(params.density*(1./PI)+.5);//poisson(u, params.density*(1./PI)); //number of impulses in nbr cell
 	vec4 h = params.sector; //annular sector
 	float a = params.a; //bandwidth
-	float aps = params.a_prime_square; //intermediate calculations for filtering
-	float filt_scale = aps*params.ainv*params.ainv;
+	float aps = im_params.a_prime_square; //intermediate calculations for filtering
+	float filt_scale = aps*im_params.ainv*im_params.ainv;
 	vec2 fpos = cpos - vec2(dnbr);//fragment position in cell space
 	
 	float acc = 0.;

			//position of impulse in cell space - uniform distribution
			vec2 ipos = vec2(nextRand(u), nextRand(u));
			//displacement to fragment
			vec2 delta = (fpos - ipos)*params.ainv;
			//impulse frequency, orientation - uniform distribution on input ranges
			float mfreq = pow(2., nextRand(u)*params.octaves);
			float ifreq = h.x*mfreq; 
			float iorientation = mix(h.z-.5*h.w, h.z+.5*h.w, nextRand(u));
			//evaluate kernel, accumulate fragment value
			vec2 mu = ifreq*vec2(cos(iorientation), sin(iorientation));
			float phi = nextRand(u); //phase - uniform dist [0, 1]
			float filt_exp = -.5*dot(mu, params.sigma_f_plus_g_inv*mu);
			acc+= filt_scale/mfreq*exp(-PI*aps*dot(delta,delta)+filt_exp)*cos(2.*PI*(dot(delta, params.filter*mu)+phi));
		}else {break;}
	}
	return acc;
 }
 
 float det2x2(mat2 m){
	return (m[0][0]*m[1][1] - m[0][1]*m[1][0]);
 }
 mat2 inv2x2(mat2 m){
	return (1./det2x2(m))*mat2(m[1][1], -m[0][1], -m[1][0], m[0][0]);
 }
 mat2 id2x2(){

//approximate poisson distribution uniform u
//from Galerne, Lagae, Lefebvre, Drettakis
const float poisson_epsilon = .001;
int poisson(inout float u, const float m){
	float u1 = nextRand(u);
	float u2 = nextRand(u);
	float x = sqrt(-2.*log(u1+poisson_epsilon))*cos(2.*PI*u2);
