C0=1; C1=1; w=0;
 function R(t)&=t^(2/(3*(1+w)));
 function C(t,dc):=1/(dc^2*R(t)^4)-1/R(t)^2  
function togalaxy(t0,t1,v,chig,win)
 global w;
 global C0;
 global C1;
 w=win;
 chi=0:0.0005:chig;
 dchi0:=v/R(t0); dchi1=v/R(t1);
 l=v*(t1-t0);
 C0=C(t0,dchi0); C1=C(t1,dchi1);
 st0=adaptiverunge("sqrt(R(y)^2+C0*R(y)^4)",chi,t0); 
 st1=adaptiverunge("sqrt(R(y)^2+C1*R(y)^4)",chi,t1);
tf0=st0[-1];
tf1=st1[-1];
vf0=R(tf0)*(R(tf0)^2+C0*R(tf0)^4)^(-1/2);
vf1=R(tf1)*(R(tf1)^2+C1*R(tf1)^4)^(-1/2);
lf=vf0*(tf1-tf0);
return [tf0,tf1,vf0,lf,lf/l]
endfunction