Этот пример показывает, как загрузить полетные данные и оценить силы G во время рейса.
Записанные данные содержат следующие параметры рейса:
угол нападения (альфа) в радианах,
угол заноса (бета) в радианах,
обозначенная скорость полета (IAS) в узлах,
тело угловые уровни (омега) в радианах/секунда,
downrange и crossrange положения в ногах, и
высота (высокий звук) в ногах.
load('astflight.mat');
Переменные MATLAB® создаются для угла нападения (альфа), угол заноса (бета), тело угловые уровни (омега) и высота (высокий звук) от записанных данных. Функция convangvel используется, чтобы преобразовать тело угловые уровни от радианов в секунду (rad/s) до степеней в секунду (градус/с).
alpha = fltdata(:,2); beta = fltdata(:,3); omega = convangvel( fltdata(:,5:7), 'rad/s', 'deg/s' ); alt = fltdata(:,10);
В этом наборе полетных данных была зарегистрирована обозначенная скорость полета (IAS). Обозначенная скорость полета (IAS) отображена в инструментировании кабины. Чтобы выполнить вычисления, истинная скорость полета (TAS), скорость полета без погрешностей измерения, обычно используется.
Погрешности измерения введены через pitot-статические анемотахометры, используемые, чтобы определить скорость полета. Эти погрешности измерения являются ошибкой плотности, ошибкой сжимаемости и калибровочной ошибкой. Применение этих ошибок к истинной скорости полета приводит к обозначенной скорости полета.
Ошибка плотности происходит должная понизить воздушную плотность на высоте. Эффект является анемотахометром чтения ниже, чем истинная скорость полета на больших высотах. Когда различие или ошибка в воздушной плотности на высоте от воздушной плотности в стандартный день на уровне моря применяются к истинной скорости полета, это приводит к эквивалентной скорости полета (EAS). Эквивалентная скорость полета является истинной скоростью полета, измененной с изменениями в атмосферной плотности, которые влияют на анемотахометр.
Ошибка сжимаемости происходит, потому что воздух имеет ограниченные возможности сопротивляться сжатию. Эта способность уменьшается увеличением высоты, увеличением скорости или ограниченным объемом. В анемотахометре существует определенное количество захваченного воздуха. При полете на больших высотах и более высоких скоростях полета, калиброванная скорость полета (CAS) всегда выше, чем эквивалентная скорость полета. Калиброванная скорость полета является эквивалентной скоростью полета, измененной с эффектами сжимаемости воздуха, которые влияют на анемотахометр.
Калибровочная ошибка характерна для данной конструкции самолета. Калибровочная ошибка является результатом положения и размещением статического вентилятора (вентиляторов), чтобы поддерживать давление, равное атмосферному давлению в анемотахометре. Положение и размещение статического вентилятора наряду с углом нападения и скоростью самолета определят давление в анемотахометре и таким образом количестве калибровочной ошибки анемотахометра. Калибровочная таблица обычно дается в экспериментальном операционном руководстве (POH) или в других спецификациях самолета. Используя эту калибровочную таблицу, обозначенная скорость полета (IAS) определяется от калиброванной скорости полета путем изменения его с калибровочной ошибкой анемотахометра.
Следующие данные являются калибровочной таблицей скорости полета для анемотахометра самолета с нулевым отклонением откидной створки. Калибровочная таблица скорости полета преобразовывает обозначенную скорость полета (IAS) в калиброванную скорость полета (CAS) путем удаления калибровочной ошибки.
flaps0IAS = 40:10:140; flaps0CAS = [43 51 59 68 77 87 98 108 118 129 140];
Обозначенная скорость полета (IAS) от рейса и калибровочной таблицы скорости полета используется, чтобы определить калиброванную скорость полета (CAS) для рейса.
CAS = interp1( flaps0IAS, flaps0CAS, fltdata(:,4) );
Атмосферные свойства, температура (T), скорость звука (a), давление (P), и плотность (ро), определяются на высоте в течение стандартного дня с помощью функции atmoscoesa.
[T, a, P, rho]= atmoscoesa( alt );
Если калиброванная скорость полета (CAS) и атмосферные свойства определяется, истинная скорость полета (Vt) может быть вычислена с помощью функции correctairspeed.
Vt = correctairspeed( CAS, a, P, 'CAS', 'TAS' );
Используйте функцию datcomimport, чтобы принести Цифровые данные DATCOM в MATLAB. Модули для этой аэродинамической информации являются ногами и степенями.
data = datcomimport( 'astflight.out', true, 0 );
Это видно в Цифровом выходном файле DATCOM и исследовании импортированных данных это
имейте данные только в первом альфа-значении. По умолчанию недостающие точки данных установлены в 99 999. Недостающие точки данных заполнены значениями для первой альфы, поскольку эти точки данных предназначаются, чтобы использоваться для всех альфа-значений.
aerotab = {'cyb' 'cnb' 'clq' 'cmq'};
for k = 1:length(aerotab)
for m = 1:data{1}.nmach
for h = 1:data{1}.nalt
data{1}.(aerotab{k})(:,m,h) = data{1}.(aerotab{k})(1,m,h);
end
end
end
Устойчивость и динамические производные в цифровой структуре DATCOM являются 3-D таблицами, которые являются функциями Числа Маха, углом нападения в градусах и высотой в ногах. В порядке выполнить 3-D линейную интерполяцию (interp3), индексы для производных таблиц требуются, чтобы быть монотонными и плед. Индексы этой формы сгенерированы функцией meshgrid.
[mnum, alp, h] = meshgrid( data{1}.mach, data{1}.alpha, data{1}.alt );
Поскольку угловые модули производных в градусах, единицы угла нападения (альфа) преобразованы от радианов от степеней функциональным convang.
alphadeg = convang( alpha, 'rad', 'deg' );
Числа Маха для рейса вычисляются функциональным machnumber с помощью скорости звука (a) и посылают (Vt) авиапочтой.
Mach = machnumber( convvel( [Vt zeros(size(Vt,1),2)], 'kts', 'm/s' ), a );
Цикличное выполнение через условия рейса, позволяет interp3 использоваться, чтобы линейно интерполировать производные таблицы, чтобы найти статические и динамические производные при тех условиях рейса.
for k = length(alt):-1:1 cd(k,:) = interp3( mnum, alp, h, data{1}.cd, Mach(k), alphadeg(k), alt(k), 'linear'); cyb(k,:) = interp3( mnum, alp, h, data{1}.cyb, Mach(k), alphadeg(k), alt(k), 'linear'); cl(k,:) = interp3( mnum, alp, h, data{1}.cl, Mach(k), alphadeg(k), alt(k), 'linear'); cyp(k,:) = interp3( mnum, alp, h, data{1}.cyp, Mach(k), alphadeg(k), alt(k), 'linear'); clad(k,:) = interp3( mnum, alp, h, data{1}.clad, Mach(k), alphadeg(k), alt(k), 'linear'); end
Если производные найдены для условий рейса, аэродинамические коэффициенты могут быть вычислены.
Ссылочные длины и области, используемые в аэродинамическом содействующем вычислении, извлечены от цифровой структуры DATCOM.
cbar = data{1}.cbar;
Sref = data{1}.sref;
bref = data{1}.blref;
Угловые модули производных в градусах, таким образом, единицы угла заноса (бета) преобразованы от радианов от степеней функциональным convang.
betadeg = convang( beta, 'rad', 'deg' );
В порядке вычислить аэродинамические коэффициенты, тело угловые уровни (омега) должны быть даны в осях устойчивости, как производные. Функциональный dcmbody2wind генерирует матрицу направляющего косинуса для осей тела к осям устойчивости (Tsb), когда угол заноса (бета) обнуляется.
Tsb = dcmbody2wind( alpha, zeros(size(alpha)) );
Скорость изменения в углу нападения (alpha_dot) также необходима, чтобы найти угловые уровни на оси устойчивости (omega_stab). Функциональный diff используется на альфе, в градусах разделенной на шаг расчета данных (0,50 секунды), чтобы аппроксимировать скорость изменения в углу нападения (alpha_dot).
alpha_dot = diff( alphadeg/0.50 );
Последнее значение alpha_dot, как сохранилось, сохраняет длину alpha_dot сопоставимой с другими массивами в этом вычислении. Это необходимо, потому что функция diff возвращает массив, который является одним значением короче что вход
alpha_dot = [alpha_dot; alpha_dot(end)];
Угловые уровни на оси устойчивости (omega_stab) вычисляются для полетных данных. Угловые уровни изменены в 3-D матрицу, которая будет умножена с 3-D матрицей для матрицы направляющего косинуса для осей тела к осям устойчивости (Tsb).
omega_temp = reshape((omega - [zeros(size(alpha)) alpha_dot zeros(size(alpha))])',3,1,length(omega)); for k = length(omega):-1:1 omega_stab(k,:) = (Tsb(:,:,k)*omega_temp(:,:,k))'; end
Вычислите коэффициент сопротивления (CD), коэффициент силы стороны (CY) и коэффициент лифта (CL). Функция convvel используется, чтобы получить модули скорости полета (Vt), сопоставимый с теми из производных.
CD = cd; CY = cyb.*betadeg + cyp.*omega_stab(:,1)*bref/2./convvel(Vt,'kts','ft/s'); CL = cl + clad.*alpha_dot*cbar/2./convvel(Vt,'kts','ft/s');
Аэродинамические коэффициенты для перетаскивания, силы стороны и лифта используются, чтобы вычислить аэродинамические силы.
Динамическое давление необходимо, чтобы вычислить аэродинамические силы. Функциональный dpressure вычисляет динамическое давление скорости полета (Vt) и плотности (ро). Функция convvel используется, чтобы получить модули скорости полета (Vt), сопоставимый с теми из плотности (ро).
qbar = dpressure( convvel( [Vt zeros(size(Vt,1),2)], 'kts', 'm/s' ), rho );
Чтобы найти силы в осях тела, матрица направляющего косинуса для осей устойчивости, чтобы придать форму оси (Tbs) необходима. Матрица направляющего косинуса для осей устойчивости, чтобы придать форму оси (Tbs) является транспонированием матрицы направляющего косинуса для осей тела к осям устойчивости (Tsb). Чтобы взять транспонирование трехмерного массива, функция permute используется.
Tbs = permute( Tsb, [2 1 3] );
Цикличное выполнение через точки полетных данных, аэродинамические силы вычислены и преобразованы от устойчивости до осей тела. Функция convpres используется, чтобы получить единицы динамического давления (qbar) сопоставимый с теми из области ссылки (Sref).
for k = length(qbar):-1:1 forces_lbs(k,:) = Tbs(:,:,k)*(convpres(qbar(k),'Pa','psf')*Sref*[-CD(k); CY(k); -CL(k)]); end
Постоянная тяга оценивается в осях тела.
thrust = ones(length(forces_lbs),1)*[200 0 0];
Постоянная оценка тяги добавляется к аэродинамическим силам, и единицы преобразованы к метрике.
forces = convforce((forces_lbs + thrust),'lbf','N');
Используя расчетные силы, оцените силы G во время рейса.
Ускорения оцениваются с помощью расчетных сил и массы, преобразованной в килограммы с помощью convmass. Ускорения преобразованы в силы G, использующие convacc.
N = convacc(( forces/convmass(84.2,'slug','kg') ),'m/s^2','G''s'); N = [N(:,1:2) -N(:,3)];
G силы построены по рейсу.
h1 = figure; plot(fltdata(:,1), N); xlabel('Time (sec)') ylabel('G Force') title('G Forces over the Flight') legend('Nx','Ny','Nz','Location','Best')
close(h1);