g_error_geofit.m 2.73 KB
``````function errGeoFit = g_error_geofit(cv,imgWidth,imgHeight,xp,yp,ic,jc,...
hfov,lambda,phi,H,theta,...
hfov0,lambda0,phi0,H0,theta0,...
hfovGuess,lambdaGuess,phiGuess,HGuess,thetaGuess,...
dhfov,dlambda,dphi,dH,dtheta,...
LON0,LAT0,...
i_gcp,j_gcp,lon_gcp,lat_gcp,...
theOrder,field);

for i=1:length(theOrder)
if theOrder(i) == 1; hfov   = cv(i); end
if theOrder(i) == 2; lambda = cv(i); end
if theOrder(i) == 3; phi    = cv(i); end
if theOrder(i) == 4; H      = cv(i); end
if theOrder(i) == 5; theta  = cv(i); end
end

% Perform the geometrical transformation
[LON LAT] = g_pix2ll(xp,yp,imgWidth,imgHeight,ic,jc,...
hfov,lambda,phi,theta,H,LON0,LAT0,field);

% Calculate the error between ground control points (GCPs)
% and image control points once georectified.
%
% This error (errGeoFit) is the error associated with the geometrical fit,
% as opposed to the error associated with a polynomial fit that can be
% calculated if requested by the user.
%

% Determine the number of CGP
ngcp = length(i_gcp);

ngcpFinite = 0;
errGeoFit  = 0.0;

for k = 1:ngcp

% Calculate the distance (i.e. error) between GCP and rectificed ICPs.
distance = g_dist(lon_gcp(k),lat_gcp(k),LON(k),LAT(k),field);

% Check if the distance is finite. The distance may be NaN for some
% GCPs that may temporarily be above the horizon. Those points are
% blanked out in the function g_pix2ll.
if isfinite(distance) == 1
errGeoFit = errGeoFit + distance^2;
ngcpFinite = ngcpFinite + 1;
else
errGeoFit = Inf;
end

end

% rms distance
errGeoFit = sqrt(errGeoFit/ngcpFinite);

% Check if the parameters are within the specified uncertainties.
% If not set the error to infinity.
if abs(hfov - hfovGuess)     > dhfov;   errGeoFit = inf; end
if abs(lambda - lambdaGuess) > dlambda; errGeoFit = inf; end
if abs(phi - phiGuess)       > dphi;    errGeoFit = inf; end
if abs(H - HGuess)           > dH;      errGeoFit = inf; end
if abs(theta - thetaGuess)   > dtheta;  errGeoFit = inf; end

%fprintf('\n')
%fprintf('Field of view (fov):             %f\n',hfov)
%fprintf('Dip angle (lambda):              %f\n',lambda)
%fprintf('Tilt angle (phi):                %f\n',phi)
%fprintf('Camera altitude (H):             %f\n',H)
%fprintf('View angle from North (theta):   %f\n',theta)
%fprintf('RMS error (m):                   %f\n',errGeoFit)
%fprintf('\n')

%figure(100)
%hold on;
%plot(errGeoFit);

``````