getuncert

 GETUNCERT finds the uncertainty in the calibration parameters

 GETUNCERT finds the uncertainty in the calibration parameters using
 jack-knife resampling method.
 
 INPUTS:
     Lpts: 3xN array containing all N points to be used for the
     calibration process.
 
     Lptsnos: 3xN array containing the corresponding scan numbers of the
     points in Lpts.
     
     Nc: 3xN array containing the corresponding normal vectors for the
     points in Lpts.
 
 INPUTS:
     deltae: 3x1 vector with the standard deviations of the delta vector.
     
     rote: 3x1 vector with the standard deviations of the Euler angles.

0001 function [deltae,rote]=getuncert(Lpts,Lptsnos,Nc)
0002 % GETUNCERT finds the uncertainty in the calibration parameters
0003 %
0004 % GETUNCERT finds the uncertainty in the calibration parameters using
0005 % jack-knife resampling method.
0006 %
0007 % INPUTS:
0008 %     Lpts: 3xN array containing all N points to be used for the
0009 %     calibration process.
0010 %
0011 %     Lptsnos: 3xN array containing the corresponding scan numbers of the
0012 %     points in Lpts.
0013 %
0014 %     Nc: 3xN array containing the corresponding normal vectors for the
0015 %     points in Lpts.
0016 %
0017 % INPUTS:
0018 %     deltae: 3x1 vector with the standard deviations of the delta vector.
0019 %
0020 %     rote: 3x1 vector with the standard deviations of the Euler angles.
0021 
0022 % setup for closed form solution
0023 Lptshat = [Lpts(1,:); Lpts(3,:); ones(1,size(Lpts,2))];
0024 
0025 Avec=zeros(size(Lpts,2),9);
0026 
0027 for cntr=1:size(Lpts,2)
0028     Avec(cntr,:)=reshape(Lptshat(:,cntr)*Nc(:,cntr)',1,9);
0029 end
0030 
0031 N=sqrt(sum(Nc.^2)');
0032 
0033 % normalise
0034 Avec=Avec./repmat(N,1,9);
0035 
0036 % get indices for each scan
0037 scannos=unique(Lptsnos);
0038 
0039 sinds=cell(1,max(scannos));
0040 
0041 for cntr=scannos
0042     sinds{cntr}=find(Lptsnos==cntr);
0043 end
0044 
0045 scancomb=combnk(scannos,length(scannos)-1);
0046 % scancomb=combnk(scannos,5);
0047 
0048 % allocate memory for parameters
0049 
0050 parmat=zeros(size(scancomb,1),6);
0051 detmat=zeros(size(scancomb,1),1);
0052 
0053 options = optimset('LargeScale','off','Display','off');
0054 warning('off','optim:fmincon:NLPAlgLargeScaleConflict');
0055 
0056 Ns=size(scancomb,1);
0057 
0058 % loop over all scan combinations
0059 for cntr=1:size(scancomb,1)
0060     ctinds=cell2mat(sinds(scancomb(cntr,:)));
0061     A=Avec(ctinds,:);
0062     ctN=N(ctinds);
0063     h=A\ctN;
0064     H=reshape(h,3,3)';
0065     ctphi=[H(:,1), cross(-H(:,1),H(:,2)), H(:,2)]';
0066     ctphi0=ctphi;
0067     ctphi = fmincon(@(ctphi)frobenius_norm(ctphi,ctphi0),ctphi0,[],[],[],[],[],[],@constraint_phi,options);
0068     ctrot=dcm2angvec(ctphi);
0069     ctdelta=H(:,3); % delta positive
0070     parmat(cntr,:)=[ctdelta',ctrot'];
0071     detmat(cntr)=1/cond(A);
0072 end
0073 
0074 % jack knifing
0075 deltamat=parmat(:,1:3);
0076 rotmat=parmat(:,4:6);
0077 
0078 % get average of delta
0079 deltamatm=sum(deltamat)/Ns;
0080 
0081 % get average of rot
0082 [rotmatx,rotmaty]=pol2cart(rotmat,ones(size(rotmat)));
0083 rotmatxm=sum(rotmatx)/Ns;
0084 rotmatym=sum(rotmaty)/Ns;
0085 rotmatm=cart2pol(rotmatxm,rotmatym);
0086 
0087 % get differences
0088 deltamatd=deltamat-repmat(deltamatm,Ns,1);
0089 rotmatd=rotmat-repmat(rotmatm,Ns,1);
0090 
0091 
0092 % adjust for angle differences more than pi and less than -pi
0093 rotmatd(rotmatd>pi)=rotmatd(rotmatd>pi)-2*pi;
0094 rotmatd(rotmatd<-pi)=rotmatd(rotmatd<-pi)+2*pi;
0095 
0096 deltae=sqrt((Ns-1)/Ns*sum(deltamatd.^2))';
0097 rote=sqrt((Ns-1)/Ns*sum(rotmatd.^2))';
0098 
0099