function twopoint;
% calculate profile by two point method with pitching calibration
% measured values are obtained from csv file
% meaning of each column must be known
% data
% 丂twop1, twop2: data of distance sensor (unit: um乯
%   pitch1: data of angle sensor乮unit: arcsec->rad乯

% parameter
%   s: sampling interval丂um
%   D: sensor interval丂um

s =1000;
D1 =25000;  
d1=D1./s;
ArcRad = pi./(3600*180);

M = dlmread('Data300.txt', '\t',0,0);
CalcuType='x mirror';

row=size(M,1);
% total data number
N_all = row;
Nd1 = N_all - d1;
r1=N_all./d1;
nn = 0:N_all-1;

if CalcuType =='y mirror'
    % 2nd column is pitching
    pitch1=M(1:Nd1,2);
    % 3rd column is distance of sensor1
    twop1=M(1:Nd1,3).*1000;
    %  4th column is distance of sensor2
    twop2=M(1:Nd1,4).*1000;
else
    pitch1=M(d1+1:N_all,2);
    twop1=M(d1+1:N_all,3).*1000;
    twop2=M(d1+1:N_all,4).*1000;
end

%display read data
figure(1);
subplot(2,1,1);
plot([1:Nd1], twop1,'-b',[1:Nd1], twop2,':r');
title('data of two distance sensors  um');
subplot(2,1,2);
plot([1:Nd1],pitch1);
title('data of angle sensor丂arcsec');

%calculate the difference of two point data
for jj=1:Nd1;
       gd1x(jj) = twop2(jj) - twop1(jj) - D1.* pitch1(jj).*ArcRad;
end;

% calculate the residual data from difference by nature extend method
% [1, Nd]
fd1x=gd1x;
% [Nd+1,N]
for ii=Nd1+1:N_all;
    fd1x(ii) = 0;
    for jj=ii-(r1-1).*d1:d1:ii-d1;
        fd1x(ii) = fd1x(ii) - gd1x(jj);
    end;
end;

% fourier transformation
fd1k = fft(fd1x);

% using transfer function H(K) to calculate f1k, f2k丂
for kk=1:N_all;
    if rem(kk-1, r1) == 0
        f1k(kk) = 0;
    else
        f1k(kk) = fd1k(kk)./(exp(2.*pi.*j.*(kk-1)./r1) -1); % f(k)=H(k)*fd(k)
    end
end;

% inverse fourier transformation
f1x=real(ifft(f1k));

% calibration profile tilt
p0 = f1x(1).*ones(1,N_all);
p1 = 0:(f1x(N_all)-f1x(1))./(N_all-1):f1x(N_all)-f1x(1);
f1x = f1x - p0 - p1;

% calculate profile by integrate method
pitch1all=M(:,2).*ArcRad;
twop1all=M(:,3).*1000;
twop2all=M(:,4).*1000;
for jj=1:N_all;
    temp = (twop2all(jj) - twop1all(jj))./d1 - pitch1all(jj).*s;
    if jj==1
        t1x(jj)=temp;
    else
        t1x(jj)=t1x(jj-1)+temp;
    end
end;
p0 = t1x(1).*ones(1,N_all);
p1 = 0:(t1x(N_all)-t1x(1))./(N_all-1):t1x(N_all)-t1x(1);
t1x = t1x - p0 - p1;

% display the profile
figure(2);
plot(nn,f1x,'-b',nn,t1x,':r');
title('profile  blue: two point; red: integrate');
mingf1=min(f1x);
maxgf1=max(f1x);
axis([-2, N_all+1, mingf1-abs(mingf1).*0.1, maxgf1+abs(maxgf1).*0.1]);

%write to file
%M=[nn;f1x;t1x];
%dlmwrite(sprintf('%s.txt',CalcuType),M','delimiter','\t','precision','%.6f');

%keyboard
return