made mod_euler work I think (not tested yet though

particle.h

@@ -194,9 +194,9 @@ namespace PARTICLE{
 	  float tmp=rand(); tmp/=RAND_MAX;
 	  sign=tmp > 0.5 ? 1:-1;
 	}
-	m_testx = m_x+sign*m_steplength/2/m_xdim*rx_new;
-	m_testy = m_y+sign*m_steplength/2/m_ydim*ry_new;
-	m_testz = m_z+sign*m_steplength/2/m_zdim*rz_new;
+	m_testx = m_x+sign*m_steplength/m_xdim*rx_new;
+	m_testy = m_y+sign*m_steplength/m_ydim*ry_new;
+	m_testz = m_z+sign*m_steplength/m_zdim*rz_new;
 	
 	if(init){
 	  m_rx_init=m_rx;

streamlines.cc

@@ -4,8 +4,22 @@
 
 namespace TRACT{
 
-
- 
+  ColumnVector mean_sph_pol(ColumnVector& A, ColumnVector& B){
+    // A and B contain th, ph f. 
+    float th,ph;
+    ColumnVector rA(3), rB(3);
+    rA<< sin(A(1))*cos(A(2)) <<sin(A(1))*sin(A(2)) <<cos(A(1));
+    rB<< sin(B(1))*cos(B(2)) <<sin(B(1))*sin(B(2)) <<cos(B(1));
+    
+    if(SP(rA,rB).AsScalar()>0)
+      cart2sph((rA+rB)/2,th,ph);
+    else
+      cart2sph((rA-rB)/2,th,ph);
+    
+    A(1)=th; A(2)=ph;
+    return A;
+  }
+  
   void read_masks(vector<string>& masks, const string& filename){
     ifstream fs(filename.c_str());
     string tmp;
@@ -197,6 +211,7 @@ namespace TRACT{
 		  ColumnVector test_th_ph_f;
 		  m_part.testjump(th_ph_f(1),th_ph_f(2));
 		  test_th_ph_f=vols.sample(m_part.testx(),m_part.testy(),m_part.testz(),m_part.rx(),m_part.ry(),m_part.rz(),pref_x,pref_y,pref_z);
+		  test_th_ph_f=mean_sph_pol(th_ph_f,test_th_ph_f);
 		  m_part.jump(test_th_ph_f(1),test_th_ph_f(2));
 		}
 	    }