用C、VB、C++或MATLAB语言编写空间后方交会程序（1）已知航摄仪内方位元素f＝153.24mm，Xo＝Yo＝0。限差0.1秒 （2）已知4对点的影像坐标和地面坐标： 影像 坐标 地面坐标 x（mm） y（mm） X(m) Y（m） Z(m) 1 -86.15 -68.99 36589.41 25273.32 2195.17 2 -53.40 82.21 37631.08 31324.51 728.69 3 -14.78 -76.63 39100.97 24934.98 2386.50 4 10.46 64.43 40426.54 30319.81 757.31

以下是一个用MATLAB语言编写的空间后方交会程序的示例代码，可以根据给定的数据进行计算和精度评定。

% 已知数据
f = 153.24; % 内方位元素
Xo = 0; % 内方位元素
Yo = 0; % 内方位元素
tolerance = 0.1; % 限差

% 影像坐标和地面坐标数据
image_coords = [-86.15 -68.99; -53.40 82.21; -14.78 -76.63; 10.46 64.43];
ground_coords = [36589.41 25273.32 2195.17; 37631.08 31324.51 728.69; 39100.97 24934.98 2386.50; 40426.54 30319.81 757.31];

% 计算外方位元素
[Xs, Ys, Zs, omega, phi, kappa] = space_resection(image_coords, ground_coords, f, Xo, Yo, tolerance);

% 输出结果
fprintf('外方位元素：\n');
fprintf('Xs = %.2f m, Ys = %.2f m, Zs = %.2f m, omega = %.2f deg, phi = %.2f deg, kappa = %.2f deg\n', Xs, Ys, Zs, rad2deg(omega), rad2deg(phi), rad2deg(kappa));

% 计算重心高程误差
e = check_point_error(image_coords, ground_coords, f, Xo, Yo, Xs, Ys, Zs, omega, phi, kappa);
fprintf('重心高程误差：%.2f m\n', e);

% 空间后方交会函数
function [Xs, Ys, Zs, omega, phi, kappa] = space_resection(image_coords, ground_coords, f, Xo, Yo, tolerance)
    % 坐标转换
    x = image_coords(:, 1) + Xo;
    y = image_coords(:, 2) + Yo;
    X = ground_coords(:, 1);
    Y = ground_coords(:, 2);
    Z = ground_coords(:, 3);
    
    % 初始值
    omega = 0;
    phi = 0;
    kappa = 0;
    Xs = mean(X);
    Ys = mean(Y);
    Zs = mean(Z);
    
    % 迭代计算
    max_iter = 100;
    iter = 0;
    while iter < max_iter
        % 计算法线方程系数矩阵
        A = [];
        for i = 1:length(x)
            row = [X(i) Y(i) Z(i) 0 0 0 -x(i)*X(i) -x(i)*Y(i) -x(i)*Z(i) x(i)];
            A = [A; row];
            row = [0 0 0 X(i) Y(i) Z(i) -y(i)*X(i) -y(i)*Y(i) -y(i)*Z(i) y(i)];
            A = [A; row];
        end
        
        % 计算误差向量
        b = [];
        for i = 1:length(x)
            b = [b; x(i) - X(i)*cos(phi)*cos(kappa) - Y(i)*(cos(omega)*sin(kappa) - sin(omega)*sin(phi)*cos(kappa)) - Z(i)*(sin(omega)*sin(kappa) + cos(omega)*sin(phi)*cos(kappa))];
            b = [b; y(i) - X(i)*cos(phi)*sin(kappa) - Y(i)*(cos(omega)*cos(kappa) + sin(omega)*sin(phi)*sin(kappa)) - Z(i)*(-sin(omega)*cos(kappa) + cos(omega)*sin(phi)*sin(kappa))];
        end
        
        % 解方程
        dx = inv(A'*A)*A'*b;
        
        % 更新外方位元素
        Xs = Xs + dx(1);
        Ys = Ys + dx(2);
        Zs = Zs + dx(3);
        omega = omega + dx(4);
        phi = phi + dx(5);
        kappa = kappa + dx(6);
        
        % 判断是否满足限差要求
        if max(abs(dx(1:3))) < tolerance && max(abs(rad2deg(dx(4:6)))) < tolerance
            break;
        end
        
        iter = iter + 1;
    end
end

% 计算重心高程误差函数
function e = check_point_error(image_coords, ground_coords, f, Xo, Yo, Xs, Ys, Zs, omega, phi, kappa)
    % 坐标转换
    x = image_coords(:, 1) + Xo;
    y = image_coords(:, 2) + Yo;
    X = ground_coords(:, 1);
    Y = ground_coords(:, 2);
    Z = ground_coords(:, 3);
    
    % 计算重心高程
    g = (Xs*cos(phi)*cos(kappa) + Ys*(cos(omega)*sin(kappa) - sin(omega)*sin(phi)*cos(kappa)) + Zs*(sin(omega)*sin(kappa) + cos(omega)*sin(phi)*cos(kappa))) / (cos(phi)*cos(kappa) + cos(omega)*sin(kappa)*sin(phi) + sin(omega)*sin(kappa)*cos(phi));
    
    % 计算重心高程误差
    e = mean(Z - g);
end

这段代码实现了空间后方交会的计算和精度评定，通过给定的影像坐标和地面坐标数据，计算出外方位元素和重心高程误差，并输出结果。