 	thrust::device_vector<double> d_theta(2);
 
 	// init, etc
-	d_theta[0] = 1;
+	d_theta[0] = 5;
 	d_theta[1] = 2;
 	// cudaMemcpyToSymbol(c_theta, p_theta, d_theta.size() * sizeof(*p_theta));
 	{



// constants for integration
static __constant__ double c_absc[] = {
	-3.4361591188377352, -2.5327316742327906, -1.756683649299881, -1.0366108297895147,
	-0.34290132722370503, 0.34290132722370431, 1.0366108297895129, 1.7566836492998816,

	if (i<n)
	{
		double mu = theta[0];
		double var = theta[1];
		double x;
		double rho;
		marg[i] = 0;
		for (int j=0; j<N_GH; j++) {
			x = features[0][i] + c_rt2*sigmas[0][i]*c_absc[j]
			rho = exp(-(x-mu)**2/(2.*var))/c_rt2pi
			marg[i] += c_wts[j]*rho;
		}
	}
}
	  
