Index: /branches/eam_branches/ipp-20151113/Ohana/src/uniphot/doc/galmotion.sh
===================================================================
--- /branches/eam_branches/ipp-20151113/Ohana/src/uniphot/doc/galmotion.sh	(revision 39205)
+++ /branches/eam_branches/ipp-20151113/Ohana/src/uniphot/doc/galmotion.sh	(revision 39206)
@@ -4,17 +4,22 @@
 
 ## FEAST-HIP
-$A_oort = +14.82; # km/sec/kpc
-$B_oort = -12.37; # km/sec/kpc
-$U_sol  =   9.32; # km/sec
-$V_sol  =  11.18; # km/sec
-$W_sol  =   7.61; # km/sec
+if (1)
+  $A_oort = +14.82; # km/sec/kpc
+  $B_oort = -12.37; # km/sec/kpc
+  $U_sol  =   9.32; # km/sec
+  $V_sol  =  11.18; # km/sec
+  $W_sol  =   7.61; # km/sec
+end
 
 ## ROESSER:
-$A_oort = +14.50; # km/sec/kpc
-$B_oort = -13.00; # km/sec/kpc
-$U_sol  =   9.44; # km/sec
-$V_sol  =  11.90; # km/sec
-$W_sol  =   7.20; # km/sec
-
+if (0)
+  $A_oort = +14.50; # km/sec/kpc
+  $B_oort = -13.00; # km/sec/kpc
+  $U_sol  =   9.44; # km/sec
+  $V_sol  =  11.90; # km/sec
+  $W_sol  =   7.20; # km/sec
+end
+
+# galactic to celestial
 $trans_phi = -62.8717488056;
 $trans_Xo  = 282.8594812080;
@@ -34,4 +39,84 @@
 # generate a collection of stars in R,D [L,B] with a range of distances
 # generate the predicted proper motions
+
+macro mkstars_uniform
+  if ($0 != 5)
+    echo "USAGE: mkstars_uniform (Nstars) (Rmax) (Dm_err) (Nsample)"
+    break
+  end
+
+  $Nstars = $1
+  $Rmax = $2
+  $dDm_o = $3
+
+  create seq 0 $Nstars
+  set phi   = 2*rnd(seq) - 1
+
+  # uniform density
+  # R ~ R_max^3/2
+  set d_kpc_o = $Rmax*rnd(seq)^(1/3)
+  set Dm = 5*log(d_kpc_o * 100)
+
+  set L = 360*rnd(seq)
+  set B = dasin(phi)
+
+  set R = L
+  set D = B
+  csystem G C R D
+
+  set uL_gal =     ($A_oort * dcos(2.0*L) + $B_oort) * dcos(B*1.0) * $iFkap;
+  set uB_gal = -0.5*$A_oort * dsin(2.0*L) *            dsin(B*2.0) * $iFkap;
+
+  set uL_sol =  ($U_sol * dsin(L) - $V_sol * dcos(L))                           * $iFkap / d_kpc_o;
+  set uB_sol = (($U_sol * dcos(L) + $V_sol * dsin(L))*dsin(B) - $W_sol*dcos(B)) * $iFkap / d_kpc_o;
+
+  set uL_o = uL_gal + uL_sol
+  set uB_o = uB_gal + uB_sol
+
+  set uL_sol_o = uL_sol
+  set uB_sol_o = uB_sol
+
+  # create a buffer to store the results : results for each star are in the y-dir, 
+  mcreate resultsL $Nstars $4
+  mcreate resultsB $Nstars $4
+
+  # now MC a number of distances
+  for i 0 $4
+    echo $i
+
+    gaussdev dDm $Nstars 0.0 $dDm_o
+    set Dm_x = Dm + dDm
+
+    set d_kpc = ten(0.2*Dm_x - 2)
+    
+    set uL_sol =  ($U_sol * dsin(L) - $V_sol * dcos(L))                           * $iFkap / d_kpc;
+    set uB_sol = (($U_sol * dcos(L) + $V_sol * dsin(L))*dsin(B) - $W_sol*dcos(B)) * $iFkap / d_kpc;
+
+    set duL = uL_sol - uL_sol_o
+    set duB = uB_sol - uB_sol_o
+
+    mset resultsL duL -x $i    
+    mset resultsB duB -x $i    
+  end
+    
+  set duL = zero(uL_o)
+  set duB = zero(uB_o)
+  set uLoff = zero(uL_o)
+  set uBoff = zero(uB_o)
+
+  for i 0 $Nstars
+    mget resultsL tmp -y $i
+    vstat -q tmp
+    
+    uLoff[$i] = $MEAN
+    duL[$i] = $SIGMA
+
+    mget resultsB tmp -y $i
+    vstat -q tmp
+    
+    uBoff[$i] = $MEAN
+    duB[$i] = $SIGMA  
+  end
+end
 
 macro mkstars_dmerr