function  [out,fT] = propsw(datin,method,Toff,startend)
%	Function to propagate solar wind at upstream monitor to another
%	spacecraft along the X axis at velocity of solar wind
%	Created: March 19, 2003
%       datin.dnup  =   datenumber at the upstream spacecraft
%       datin.Rup   =   the GSM coordinates of the upstream spacecraft in Re [X, Y, Z]
%       datin.Bup   =   the GSM components  of the upstream magnetic field
%       datin.Rdn   =   the GSM coordinates of the downstream spacecraft
%       datin.Vup   =   the solar wind speed at the upstream spacecraft [|V|
%       datin.dat   =   a column array for which every column will be delayed (all sampled at times in dnup)
%       out         =   [dnnew, propagated data]    This is upstream time plus delay and prop data interp to new grid        

dn      = datin.dnup;                               %Pick out the time at upstream spacecraft
dXgsm   = datin.Rup(:,1)-datin.Rdn(:,1);            %Total X distance in Re
dYgsm   = datin.Rup(:,2)-datin.Rdn(:,2);            %Total Y distance in Re
if(strcmp(method,'Spherical'))
    tanpsi  = 0;                                        %Assume Fronts orthog to radius
    T = 6378*(dXgsm+tanpsi*dYgsm)./datin.Vup;         %Time delay in seconds
elseif(strcmp(method,'Parker'))
    tanpsi = 1;
    T = 6378*(dXgsm+tanpsi*dYgsm)./datin.Vup;         %Time delay in seconds
elseif(strcmp(method,'Spiral'))
    alp = Rhodes04(datin.Bup);     
    tanpsi = tan(alp);
    T = 6378*(dXgsm+tanpsi.*dYgsm)./datin.Vup;        %Time delay in seconds
elseif(strcmp(method,'Block'))
    T   = repmat(Toff,length(dn),1);                  %Make an array of fixed offset
else
    disp('No valid method specified')
end
fT      = ConvMean(T/60,-7,7,.1);                       %Filtered time delay in minutes
fT(1:7) = T(1:7)/60;                                    %Replace the nan
fT(end-6:end) = T(end-6:end)/60;
fT      = fT/1440;                                      %Need it as fraction of a day
%There are still nan left in my test set I remove them below

% ADD DELAY TIME TO EVERY WIND TIME
tplus	= repmat(nan,length(dn),1);               %Array to hold the delayed time
is		= find(~isnan(fT));							    %Find any remaining flags in low pass tdel
tplus(is) = dn(is)+fT(is);							    %Estimated arrival time at Wind

% These are the variables to propagate
dat = datin.dat;
[md,nd] = size(dat);

% Remove all records for which tplus is nan
ij      = find(isnan(tplus));
tplus(ij) = [];
dat(ij,:)   = [];

% CREATE THE FIRST DIFFERENCE SERIES OF TPLUS AND LOCATE AND REMOVE ALL NEGATIVE TIME STEPS
dtplus	= diff(tplus);
is		= find(dtplus<=0);								    %Get all negative time changes

% Remove any points overtaken by faster moving plasma
while(length(is)>0)
   tplus(is+1) 	= [];									    %Throw out overtaken time, pressure, vel
   dat(is+1,:)  = [];                           
   dtplus		= diff(tplus);
   is			= find(dtplus<=0);
   fprintf('#overlaps = %d\n',length(is))
end

% INTERPOLATE PROPAGATED VELOCITY TO 1-MINUTE GRID
delt    = 1/1440;                                           %I have fixed the sample rate at 1 min
tmin    = nanmin(tplus);
t0      = floor(nanmin(tplus));                             %Start of this day
tb      = delt*(fix((tmin-t0)/delt)+1)+t0;                  %This is first sample later than start of prop data
tmax    = nanmax(tplus);
te	    = delt*(fix((tmax-t0)/delt))+t0;			        %This is first sample before last point in prop solar wind
ti	    = (tb:delt:te)';									%Times at each minute
%keyboard
% Degap and interpolate the propagated solar wind variables
for j = 1:nd
    out     = degapint([tplus,dat],nan,60/86400,startend,5/1440);
end
return