%DepthData=[];
%d=((((AllData(:,9))*2-0.5)*70/4)-1)*10;
%DepthData=[DepthData;AllData(:,1) d];
close all
clear all
load DepthData   

DepthData(isnan(DepthData(:,2)),:)=[];
DepthData(isnan(DepthData(:,1)),:)=[];

figure;
plot(DepthData(:,1),DepthData(:,2))
datetick('x','HH:MM','keepticks','keeplimits')
title('Raw BED Depth Data');
xlabel('Time');
ylabel('Depth (m)');

%ThalDepth=data(:,3).*-1;
load Thalwag;
figure; plot(ThalDepth);
title('Thal Depth Records');
ylabel('Depth (m)');


%ThalDist=data(:,6);
%figure;plot(ThalDist);
close all;
Thaltime= ThalDist/3.3;     %3.3m/sec is average speed

EpochSecs=(DepthData(:,1)-719529)*86400;
BedsElapsedTime=zeros(size(EpochSecs));
for i=1:size(EpochSecs,1),
        BedsElapsedTime(i,1)=(EpochSecs(i,1)-EpochSecs(1,1))/60/60;
end

figure; plot(Thaltime/60, ThalDepth)
hold on;
plot(BedsElapsedTime.*60,DepthData(:,2),'r')
plot(Thaltime/60,ThalDepth+7.6845,'g');
legend('Thalwag','BED','Thal+offset','Location','Best');xlabel('Elapsed Min');ylabel('Depth (dBar)');


[B,I,J] = unique(round(DepthData(:,2)),'stable');
BedsRounded=[BedsElapsedTime(I),B];

figure; plot(BedsElapsedTime.*60,DepthData(:,2),'b');hold on;
plot(BedsRounded(:,1).*60,BedsRounded(:,2),'g')
xlabel('Elapsed Time(min)');ylabel('Depth (dbar)');
legend('Raw Data', 'Rounded Data','Location','Best');

[C,H,K]=unique(round(ThalDepth),'stable');

figure; plot(Thaltime/60, ThalDepth)
hold on;
plot(Thaltime(H)/60,C,'r--');
xlabel('Elapsed Time(min)');ylabel('Depth (dbar)');
legend('Thal depth',' Rounded Thal depth','Location','Best');
ThalRound=[C, ThalDist(H)];


[C,IA,IB] = intersect(BedsRounded(:,2),ThalRound(:,1),'rows');

BedsDist=ones(size(BedsRounded,1),1);
BedsDist(IA)=ThalRound(IB,2);
[i,j]=find(BedsDist==1);
BedsDist(i)=nan;

BedsRounded(isnan(BedsDist),:)=[];
BedsDist(isnan(BedsDist))=[];

diffDist=diff(BedsDist);
diffTime=diff(BedsRounded(:,1)*3600);
BedsTime=BedsRounded(:,1)*3600;

[i,j]=find(diffTime<0);
diffDist(i)=[];
diffTime(i,:)=[];
BedsTime(i,:)=[];

[i,j]=find(diffDist<0);
diffDist(i)=[];
diffTime(i,:)=[];
BedsTime(i,:)=[];

figure;subplot(2,1,1);bar(diffDist);ylabel('Delta Dist (m)');
subplot(2,1,2);bar(diffTime);ylabel('Delta Time (sec)');

deltas=[diffTime diffDist BedsTime(1:end-1)];
deltas(isnan(deltas(:,1)),:)=[];
deltas(isnan(deltas(:,2)),:)=[];
figure;subplot(2,1,1);bar(deltas(:,2));ylabel('Delta Dist (m)');
subplot(2,1,2);bar(deltas(:,1));ylabel('Delta Time (sec)');



EstVelocity=deltas(:,2)./deltas(:,1);


figure; bar(EstVelocity,'b');
ylabel('Velocity (m/sec)');
%figure;plot(deltas(:,3),EstVelocity);


figure;subplot(2,1,1);
plot(diff(BedsRounded(:,1)).*3600);
ylabel('seconds')
subplot(2,1,2);
plot(diff(BedsDist));
ylabel('meters')


%%%% try griding data and filtering
% X = 0:10; V = sin(X); Xq = 0:.25:10;
% Vq = interp1(X,V,Xq); plot(X,V,'o',Xq,Vq)

X=BedsElapsedTime.*3600;    %time in elapsed sec
V=DepthData(:,2);
Xq=linspace(0,X(end),size(X,1));

Vq=interp1(X,V,Xq);
figure;plot(X,V,'o',Xq,Vq)

Vfilt=medfilt1(Vq,10);

figure;
plot(X/60,V,'b',Xq/60,Vq,'g',Xq/60,Vfilt,'r'); hold on; 
legend('Raw Data','Regrid Data','Filtered Data','Location','Best');
xlabel('Elapsed Min');
ylabel('Depth (dBar)');

figure;
plot(Thaltime/60, ThalDepth);hold on; plot(Xq/60,Vfilt,'r');
legend('Thal Const vel','Filtered Data','Location','Best');
xlabel('Elapsed Min');
ylabel('Depth (dBar)');

Vround=round(Vfilt)';
Tround=round(ThalDepth);
[C,IA,IB] = intersect(Vround,Tround,'rows');

figure;plot(Vround(IA),'o');hold on;plot(Tround(IB),'xr');
legend('Beds','Thalwag','Location','Best');
title('Intersect points');
xlabel('Data point');
ylabel('Depth (dBar)');

VelEst=[];
VelEst(:,1)=Xq(IA);
VelEst(:,2)=Vround(IA);
VelEst(:,3)=ThalDist(IB);

figure;plot(VelEst(:,1)/60,VelEst(:,2));hold on; plot(X/60,V,'r');
legend('Velocity Data Set','Raw Data','Location','Best');
ylabel('Depth (dBar)');
xlabel('Elapsed Min');

figure;subplot(2,1,1);
bar(diff(VelEst(:,1)));
ylabel('seconds');
title('Delta time');

subplot(2,1,2);
bar(diff(VelEst(:,3)));
ylabel('meters')
title('Delta Dist');

a=diff(VelEst(:,3));
b=diff(VelEst(:,1));
Velocity2=diff(VelEst(:,3))./diff(VelEst(:,1));


figure; bar(Velocity2,'b');
ylabel('Velocity (m/sec)');
hold on;
bar(EstVelocity,'r');
legend('Filtered Est', 'Rounded Est','Location','Best');

%BedsRounded=[BedsElapsedTime(I),B];

[D,L,M]=unique(Vround,'stable');
UniqFiltBeds=[Xq(L)',D];
figure; plot(X,V,'b');hold on;
plot(UniqFiltBeds(:,1),UniqFiltBeds(:,2),'gx')
xlabel('Elapsed Time(min)');ylabel('Depth (dbar)');
legend('Raw Data', 'Rounded Data','Location','Best');

[C,IA,IB] = intersect(UniqFiltBeds(:,2),ThalRound(:,1),'rows');

figure;plot(UniqFiltBeds(IA,2),'o');hold on;plot(ThalRound(IB),'xr');
legend('Beds','Thalwag','Location','Best');
title('Intersect points');
xlabel('Data point');
ylabel('Depth (dBar)');

Vel3=[];
Vel3(:,1)=UniqFiltBeds(IA,1);
Vel3(:,2)=UniqFiltBeds(IA,2);
Vel3(:,3)=ThalRound(IB,2);

figure; plot(X/60,V);hold on;plot(Vel3(:,1)/60,Vel3(:,2),'r');
xlabel('Elapsed Time(min)');ylabel('Depth (dbar)');
legend('Raw Data', 'Matched Data','Location','Best');

VelEst3=diff(Vel3(:,3))./diff(Vel3(:,1));

figure; bar(VelEst3);
[i,j]=find(VelEst3>6);
VelWhacked=VelEst3;
TimeWhacked=Vel3(2:end,1);
TimeWhacked(i)=[];
VelWhacked(i)=[];
figure;bar(VelWhacked);
figure;plot(TimeWhacked,VelWhacked);

figure;subplot(2,1,1);
bar(diff(Vel3(:,1)));
ylabel('seconds');
title('Delta time');

subplot(2,1,2);
bar(diff(Vel3(:,3)));
ylabel('meters')
title('Delta Dist');

figure; [AX,H1,H2]=plotyy(X,V,Vel3(2:end,1),VelEst3);

