pindx = [1:3 7:10];
tindx = [4:6];

pointfilenames = ['image_points_02.txt';
                  'image_points_03.txt';
                  'image_points_04.txt'];
siz = size(pointfilenames);
imgnum = siz(1);

B = load('world_points.txt');
X_ext = B(pindx,:)';
wp = B(tindx,:)';

filenames = {'intrinsic_parallel.txt',
             'intrinsic_planes.txt', %principal point is modified
                                     %for zero-origin coordinate
             'intrinsic_indoor.txt',
             'intrinsic_outdoor.txt'};
siz = size(filenames);
knum = siz(1);

for k = 1:knum
  fprintf(1, 'K: %s\n', filenames{k});
  K = load(filenames{k});
  fc = [K(1,1); K(2,2)];
  cc = [K(1,3); K(2,3)];
  kc = [0;0;0;0;0];
  alpha_c = 0;

  rms = 0;
  rmsf = 0;
  for m = 1:imgnum
    A = load(pointfilenames(m,:));
    x_ext = A(pindx,:)';
    pr = A(tindx,:)';

    % Estimate extrinsic parameters
    extrinsic_computation2;
    rmsf = rmsf + sum(sum(err_reproj .* err_reproj));
    
    [ip] = project_points2(wp, omc_ext, Tc_ext, fc, cc, kc, alpha_c);
    dv = ip - pr;
    rms = rms + sum(sum(dv .* dv));
  end
  rmsf = sqrt(rmsf / imgnum / length(pindx));
  fprintf(1, 'Fitting Error: RMS = %3.6f\n', rmsf);
  rms = sqrt(rms / imgnum / length(tindx));
  fprintf(1, 'Reprojection Error: RMS = %3.6f\n', rms);
end

