Методы цифрового моделирования

Методы цифрового моделирования

ЗАДАНИЕ

На выполнение курсовой работы

«Математические методы и модели конструкторско-технологического проектирования»

. Выполнить моделирование сигнала s(t) = (1 + a∙cos2πΩt)∙(sinπf0 + sin2πf0), построить s(t), спектр сигнала, где: а = 0.75, Ω = 20 Гц, f0 = 150 Гц.

. Построить математическую модель системы (см. рис. 2.1), рассчитать и построить графики амплитудно-частотных характеристик каждого из блоков и всего устройства в целом. ЧХ блока «Т» задана таблично (см табл. 1). Выполнить моделирование работы системы при входном сигнале s(t). Построить графики выходного сигнала и его спектра.

Таблица 1.1 - АЧХ и ФЧХ блока «Т»

. Выполнить моделирование шума n(t) с экспоненциальным законом распределения вероятностей с математическим ожиданием mn = 0.3. Оценить среднее значение и дисперсию отсчетов n(t), построить гистограмму и проверить адекватность модели n(t) по критерию Пирсона.

. Выполнить моделирование работы системы при входном сигнале n(t). Построить графики спектра, автокорреляционной функции и гистограммы выходного сигнала. Выбрать статистическую модель для выходного сигнала, найти оценки ее параметров.

. Оформить расчетно-пояснительную записку согласно ДСТУ3008-95.

математический цифровой моделирование

ВВЕДЕНИЕ

Цель работы: изучить различные методы цифрового моделирования, в том числе статистического; получить практические навыки реализации алгоритмов цифрового моделирования; закрепить знания о методах обработки статистических данных.

Объекты исследования: детерминированный сигнал; математическая модель системы; шум, распределенный по экспоненциальному закону; работа системы при воздействии на нее шума.

Метод исследования - статистическое моделирование в системе Fortran.

Ожидаемые результаты: моделирование всех объектов исследования; получение графиков сигналов и процессов; выбор статистической модели системы.

1 МОДЕЛИРОВАНИЕ СИГНАЛА НА ВХОДЕ

.1 Расчет параметров моделирования

Задан сигнал s(t) = (1 + a ∙ cos2πΩt) ∙ (sinπf0 + sin2πf0). Раскрыв скобки и применив тригонометрическую формулу: sinA ∙ cosB = sin(A - B) + sin(A + B), получаем

s(t) = sin2(πf0/2)t + sin2πf0t +  sin2π()t +  sin2π()t +

+ sin2π(f0 - Ω)t + sin2π(f0 + Ω)t. (1.1)

Подставим численные значения параметров а, Ω, f0 и получим

s(t) = sin2π75t + sin2π150t + 0,375sin2π55t + 0,375sin2π95t + 0,375sin2π130t +

 +0,375sin2π170t(1.2)

Было выбрано наибольшее значение частоты - 170 Гц. Теперь используется теорема Котельникова, которая гласит, что, если аналоговый сигнал x(t) имеет ограниченный спектр, то он может быть восстановлен однозначно и без потерь по своим дискретным отсчётам, взятым с частотой более удвоенной максимальной частоты спектра Fmax:

,(1.3)

где Fmax - верхняя частота в спектре, или (формулируя по-другому) по отсчётам, взятым с периодом чаще полупериода максимальной частоты спектра;

k - коэффициент запаса и k = 20.

Время дискретизации ∆t рассчитывается по формуле:

∆t=1/ Fd.(1.2)

После расчета получаем ∆t = 2,9412 ∙ 10-4 с.

Время наблюдения процесса рассчитывается так:

,(1.3)

Количество отсчетов N (объем выборки) связано с необходимым временем наблюдения T за процессом:

N=T/∆t.(1.4)

Для анализа спектрального состава процессов используют, как правило, процедуру быстрого преобразования Фурье (БПФ), которая в Fortran обозначается как FFT(*). Существует прямое и обратное преобразование Фурье. С помощью прямого получают спектр, а с помощью обратного сигнал восстанавливают из спектра. Особенностью процедуры БПФ является требование к объему выборки N= 2м, где М - целое число. Поэтому количество отсчетов N берем равное 1024. После этого необходимо пересчитать время дискретизации по времени ∆t, которое равно 1,9531 ∙ 10-4 с.

1.2 Описание алгоритма

Для выделения мнимой составляющей в Фортране используется функция aimag, для выделения реальной - функция real. Амплитудный спектр получают с использованием cabs - взятие модуля из комплексного числа. Фазовый спектр получают с заданной точностью 10-3, для его определения используется функция atan2, которая берется от комплексного массива сигнала. Для построения всех спектров используется М = N/2 точек, т. к. спектры обладают свойством зеркальности.

В результате моделирования получены следующие графики:

а) сигнал на входе системы;

б) амплитудный спектр входного сигнала;

в) фазовый спектр входного сигнала;

г) мнимая составляющая входного сигнала;

д) реальная составляющая входного сигнала.

Рисунок 1.1 - Моделируемый входной сигнал

Рисунок 2.2 - Амплитудно-частотная характеристика входного сигнала

Рисунок 2.3 - Фазо-частотная характеристика входного сигнала

Рисунок 2.4 - Мнимая составляющая спектра входного сигнала

Рисунок 2.5 - Реальная составляющая спектра входного сигнала

2 ПОСТРОЕНИЕ МАТЕМАТИЧЕСКОЙ МОДЕЛИ СИСТЕМЫ. МОДЕЛИРОВАНИЕ ЕЕ РАБОТЫ ПРИ ВХОДНОМ СИГНАЛЕ

.1 Описание математической модели системы и ее составляющих

На рисунке 2.1 изображена блок-схема заданной системы:

Рисунок 2.1 - Блок-схема системы

Составляющие систему блоки имеют следующий вид:

а) - схема блока 9

б) - схема блока 11

в) - схема блока 18

Рисунок 2.2 - Блоки 9, 11, 18

Математическое описание блоков - это их передаточные функции:

а) для блока 9

,(2.1)

б) для блока 11

,  (2.2)

в) для блока 18

,(2.3)

где  = j∙ω - оператор Лапласа;

L - индуктивность;

R - сопротивление;

С - емкость.

Для блока Т АЧХ и ФЧХ заданы таблично (см. табл. 1.1). Для построения этих характеристик используется полином Лагранжа

,(2.4)

где хi - значения в узлах интерполяции.

Для получения графиков АЧХ и ФЧХ Т-блока используется подпрограмма интерполяции.

Общая передаточная функция была рассчитана в несколько этапов с применением правил для соединений блоков:

а) блок 11 и Т-блок соединены последовательно - их передаточные характеристики перемножаются;

W1 = W11∙WT,(2.5)

б) новый полученный блок № 1 и блок № 18 - это цепь с отрицательной обратной связью;

W2 = (2.6)

в) новый полученный блок № 2 и блок № 9 соединены последовательно - их передаточные характеристики перемножаются.

Wобщ = W2∙W9,(2.7)

где Wi - передаточная функция.

2.2 Алгоритм моделирования

Для того, чтобы получить спектр выходного сигнала нужно

,(2.4)

где W(p) - передаточная функция;

X(p) - спектр воздействия, найден с помощью прямого преобразования Фурье.

Для того, чтобы получить выходной сигнал, используется обратное преобразование Фурье. Сигнал получается на М точек, а нужно на N, поэтому достраивается “зеркальная часть” передаточной характеристики с использованием функции conjg - комплексное сопряжение.

В результате моделирования были получены следующие графики:

а) передаточные характеристики всех блоков (9, 11, 18) (см. рис. 2.6 - 2.8);

б) интерполированные АЧХ и ФЧХ Т-блока (см. рис. 2.9 - 2.10);

в) общая передаточная характеристика (см. рис. 2.11);

г) сигнал на выходе системы (см. рис. 2.12);

д) спектры выходного сигнала - амплитудный, фазовый, мнимой и реальной составляющей (см. рис. 2.13 - 2.16).

Рисунок 2.6 - Передаточная характеристика блока 9

Рисунок 2.8 - Передаточная характеристика блока 18

Рисунок 2.9 - Интерполированная АЧХ Т-блока

Рисунок 2.9 - Интерполированная ФЧХ Т-блока

Рисунок 2.10 - Общая передаточная характеристика

Рисунок 2.11 - Сигнал на выходе системы

Рисунок 2.12 - Амплитудный спектр выходного сигнала

Рисунок 2.13 - Фазовый спектр выходного сигнала

Рисунок 2.14 - Мнимая составляющая выходного сигнала

Рисунок 2.15 - Реальная составляющая выходного сигнала

3 МОДЕЛИРОВАНИЯ ШУМА С ЗАДАНОЙ ПЛОТНОСТЬЮ РАСПРЕДЕЛЕНИЯ ВЕРОЯТНОСТЕЙ

Заданный шум распределен по экспоненциальному закону.

n(t) = λ∙e-λ∙t,(3.1)

где λ - параметр.

Для моделирования задано математическое ожидание МХ = 0,3. Параметр λ связан с математическим ожиданием

МХ = 1/ λ.(3.2)

Из (3.2) был найден параметр λ = 3,333.

Для имитации реализаций непрерывных случайных величин применяют метод обратной функции, позволяющий получить из реализаций ri - стандартной случайной величины реализации yi случайной величины с требуемым законом распределения. Для вычисления стандартной случайной величины был применен мультипликативный алгоритм.

После преобразований получаем генератор непрерывных случайных величин, распределенных по экспоненциальному закону. Обратная функция имеет вид

Xi = - .(3.2)

Шум имеет вид, показанный на рисунке 3.1.

Рисунок 3.1 - Шум, распределенный по экспоненциальному закону

Для полученных непрерывных случайных величин рассчитана оценка дисперсии

,(3.3)

оценка среднеквадратического отклонения

,(3.4)

Результаты расчета:  = 9,2 ∙ 10-2, σ = 0,3.

Для построения гистограммы (см. рис. 3.2) используется подпрограмма Rhist.

Выбрано для построения 15 интервалов.

Рисунок 3.2 - Гистограмма входного шума

Для проверки согласия теоретического и экспериментального распределений используется критерий Пирсона. По значению χ2 и количеству степеней свободы τ = 13 определяется вероятность того, на сколько согласуются вышеуказанные распределения.

,(3.5)

где NPk - теоретическое количество попаданий в интервал гистограммы;

Gk - эксперементально полученое количество попаданий в интервал гистограммы.

Для того, чтобы получить теоретическое количество попаданий в n-й интервал гистограммы, были посчитаны площади каждого интервала с использованием формулы Симпсона. Далее площадь каждого интервала была умножена на количество отсчетов (N = 1024).

После расчета получаем χ2 = 11,68 и вероятность соответствия теоретического закона распределения и экспериментально полученного 0,5.

На рисунке 3.3 изображен график экспоненты и гистограмма. Видно, что они почти согласуются. Это же подтверждается расчетом χ2.

Рисунок 3.3 - Сравнение теоретического закона распределения и экспериментально полученного

4 ВЫПОЛНЕНИЕ МОДЕЛИРОВАНИЯ РАБОТЫ СИСТЕМЫ ПРИ ВХОДНОМ СИГНАЛЕ N(t)

Для того чтобы получить выходной спектр нужно передаточную функцию системы Wобщ(p) (найденную в пункте 2) умножить на спектр входного сигнала (шума) N(p). Спектр сигнала n(t), находится с помощью прямого преобразования Фурье.

,(4.1)

Далее находятся спектры: амплитудный, фазовый, мнимой и реальной составляющих (см. рис. 4.1 - 4.4).

Рисунок 4.1 - Амплитудный спектр выходного шума

Рисунок 4.2 - Фазовый спектр выходного шума

Рисунок 4.3 - Реальная составляющая выходного шума

Рисунок 4.4 - Мнимая составляющая выходного шума

С помощью быстрого обратного преобразования Фурье от спектра выходного сигнала Y(p) находим сам выходной сигнал (см. рис.4.5).

Рисунок 4.5 - Шум на выходе системы

Для того, чтобы получить автокорреляционную функцию системы (АКФ), используется теорема Винера-Хинчина. Эта теорема связывает между собой квадрат спектральной плотности сигнала и его корреляционную функцию через преобразование Фурье. По спектру выходного шума вычисляется энергетический спектр (формально это умножение на комплексно сопряженный). А после применения обратного преобразования Фурье получаем АКФ (см. рис. 4.6). Для построения ее графика было взято N/4 отсчетов.

Рисунок 4.6 - Автокорреляционная функция

Построение гистограммы выходного шума (см. рис. 4.7) проводится с помощью подпрограммы Rhist. Взято, как и в предыдущем задании 15 интервалов, для которых рассчитывается экспериментальное попадание в интервал с использованием подпрограммы для подсчета их площади по формуле Симпсона.

Рисунок 4.7 - Гистограмма выходного шума

,(4.2)

где ni - отчеты выходного сигнала;

mn - оценка математического ожидания;

σ - среднеквадратическое отклонение

,(4.3)

где DX - дисперсия.

Эксцесс показывает, на сколько островершинный полученный закон распределения и рассчитывается как

.(4.4)

Результаты расчета: mn=-2,52·10-5, σ = 3,64 · 10-2, ά = 0,32, ε = 2,32. После нанесения этих точек на плоскость моментов (см. приложение А), видно, что ближайшими к полученному закону распределения являются распределения Релея и нормальный закон.

Плотность распределения вероятностей нормального закона имеет следующий вид

,(4.5)

Закон распределения Релея имеет следующий вид

,(4.6)

где х > 0.

Расчет теоретического количества попаданий отчетов случайного процесса в интервал гистограммы проводилось с помощью подпрограммы подсчета площади по формуле Симпсона от плотности распределения нормального закона распределения и распределения Рэлея.

Однако после проверки этих законов на критерий Пирсона, получаем:

а) для распределения Рэлея χ2 = 3601;

б) для нормального распределения χ2 = 1385.

Эти значения дают вероятность соответствия меньше 0,05.

На рис. 4.8 изображено совмещение теоретического количества попаданий в интервал гистограммы для распределения Рэлея и гистограммы выходного шума, а на рис. 4.9 совмещено теоретического количества попаданий в интервал гистограммы для нормального закона распределения и гистограммы выходного шума

Рисунок 4.8 - Совмещение теоретического количества попаданий в интервал гистограммы для распределения Рэлея и гистограммы выходного шума

Рисунок 4.9 - Совмещение теоретического количества попаданий в интервал гистограммы для нормального распределения и гистограммы выходного шума

ВЫВОДЫ

Заданный входной сигнал состоит из шести синусоид. В данном сигнале имеются по шесть мнимых и реальных отрицательных составляющих на частотах 75 Гц, 150 Гц, 55 Гц, 95 Гц, 130 Гц, 170 Гц. Однако реальных составляющих быть не должно. Их появление объясняется погрешностью вычисления. К тому же их значения достаточно малы - наибольшее из них -1,963 ∙10-3. В амплитудном спектре шесть составляющих, как это и должно быть, а их амплитуды соответствуют коэффициентам при синусоидах в сигнале. В фазовом спектре все составляющие компоненты с одинаковой фазой π/2.

В связи с тем, что передаточная функция 18 блока имеет резонанс, общая передаточная характеристика имеет тоже резонансный тип. Этот же фактор влияет на вид выходного сигнала - появляются резонансные периодические выбросы. Амплитудный спектр имеет шесть составляющих. В фазовом и мнимом спектре также по шесть составляющих и все они имеют отрицательные значения. В реальном спектре появляются три отрицательные и три положительные составляющие.

Заданный шум распределен по экспоненциальному закону. При построении совместного графика теоретического распределения и гистограммы видно, что они хорошо согласуются. Это же подтверждается расчетами: χ2 = 11,68 и вероятность соответствия теоретического закона распределения и экспериментально полученного 0,5.

При прохождении шума через систему его форма искажается: форма его распределения уже не похожа на экспоненциальную. Для того, чтобы определить на какой закон распределения похож полученный на выходе шум, рассчитываются коэффициенты эксцесса и асимметрии. Далее эти точки наносятся на плоскость моментов, и определяется, какой из законов распределения ближайший. Наиболее подходящим оказался закон распределения Релея, но при проверке данного распределения по критерию Пирсона оказалось, что вероятность соответствия между экспериментально полученным законом распределения и распределением Релея меньше 0,05. это связано с тем, что Распределение Релея начинается с нуля, в отличии от экспериментально полученного. Поэтому также выполнилась оценка по критерию Пирсона на соответствие с нормальным законом распределения. В этом случае χ2 = 1385, что тоже дает вероятность соответствия меньше 0,05.

ПРИЛОЖЕНИЕ А

ПРИЛОЖЕНИЕ Б

!ПЕРВОЕ ЗАДАНИЕ - МОДЕЛИРОВАНИЕ СИГНАЛА S(t)=(1+а*соs2*pi*w*t)*(sinpi*f0*t+sin2*pi*f0*t),

!ПОСТРОЕНИЕ ГРАФИКА S(t) И СПЕКТРОВ СИГНАЛОВ.

INTEGER,PARAMETER:: N=1024, M=512, k=11dt/1.9531e-4/,df/5.0/S(N), t(N)SPECTR_A(M), SPECTR_Ph(M), SPECTR_Re(M), SPECTR_Im(M),f(M)cS(N)a/0.75/, f0/150.0/, w/20.0/pi/3.1415/R1/5.1/, R2/0.1/, C1/68.0e-6/, L1/1.0e-2/R_1/10.0/, R_2/2.0/, L_1/5.0e-2/R__1/1000.0/, R__2/1.0/, C__1/47.0e-6/,L__1/1.0e-2/W9(M),W11(M),W18(M)FT(K)/0.0,260.0,520.0,780.0, 1040.0,1300.0,1560.0,1820.0,2080.0,2340.0,2600.0/AChH(K)/0.0, -6.54, -11.5, -15.8, -18.3, -19.9, -21.7, -24.0, -25.2, -28.2, -32.1/FChH(K)/0.0, -58.0, -72.7, -78.3, -81.2, -83.0, -84.2, -85.0, -85.7, -86.2, -86.6/AChHT(M), FChHT(M), AChHT_EDENICY(M), FChHT_RADIANY(M), WTT(M), W_USTROJSTVA(M)SpA2(M),SpPh2(M),SpRe2(M),SpIm2(M) ,SY(N), PARAMETER::Nhist=15,NPOINT=50*NHIST+1,Ng=2*4 R, XSTART/0.37290981/,XKOEF/37.0/, E/2.7/,X(N),H(NHIST)MX/0.3/, DX, DX0, SIGMAX, XMAX, XMIN, TT(NHIST),FX(NPOINT),LYAMDA/3.3333/,XX(NPOINT)PLOSHA(NHIST),XF(120),XXX(NHIST),PLOSHA1(NHIST),XXXX(120),INTERVAL,DVA(NHIST,NG),TT1(NHIST,NG),XIXIXMAX2, XMIN2,TT2(NHIST),SxvyhR(N),RAKF(N),C,RAKF1(N),MXSh,DXSh, ALFAX, EX,SIGMAXSh,MX0Sh,DX0Sh,ALFAX0,EX0,H1(NHIST)XNORM(NPOINT), BNORM,ANORM,PLOTN(NPOINT),INTERVAL2,POPADANIE1(NHIST),DVA2(NHIST,NG),POPADANIE(NHIST),XNORMAL(NHIST,NG)XREL(NPOINT), BREL,AREL,PLOTN2(NPOINT),POPADANIE2(NHIST),DVA3(NHIST,NG),POPADANIE3(NHIST),XRELEJ(NHIST,NG),INTERVAL3CP(N), CW9(N),CW11(N),CW18(N), CWT(M)CW1(M), CW2(M), CW3(N), CW(N), CSY(N),CX(N),Sxvyh(N),SPECTRshuma_Re(M),SPECTRshuma_Im(M),SPECTRshuma_A(M),SPECTRshuma_Ph(M)CJ/(0.0,1.0)/,AKF(N)I=1,N(I)=REAL(I-1)*dt

S(I)=(1+a*COS(2*pi*w*t(I)))*(SIN(pi*f0*t(I))+SIN(2*pi*f0*t(I)))DO

!CALL DRAW(N,1,t,S,0,'LINE','МОДЕЛИРУЕМЫЙ ВХОДНОЙ СИГНАЛ')

!ACCEPT*I=1,N(I)=CMPLX(S(I),0.0)DOFFT (cS,N,0)I=1,M_Re(I)=REAL(cS(I))

SPECTR_Im(I)=AIMAG(cS(I))

f(I)=REAL(I-1)*dfDO

!CALL DRAW (M,1,f,SPECTR_Re,0,'STICK','РЕАЛЬНАЯ СОСТАВЛЯЮЩАЯ ВХОДНОГО СИГНАЛА')

!ACCEPT*

!CALL DRAW (M,1,f,SPECTR_Im,0,'STICK','МНИМАЯ СОСТАВЛЯЮЩАЯ ВХОДНОГО СИГНАЛА')

!ACCEPT*I=1,M_A(I)=CABS(cS(I))

IF(ABS(REAL(cS(I)))<1.0e-3.AND.ABS(AIMAG(cS(I)))<1.0e-3) THEN

SPECTR_Ph(I)=0.0

ELSE_Ph(I)=ATAN2(AIMAG(cS(I)),REAL(CS(I)))DO

!CALL DRAW (M,1,f,SPECTR_A,0,'STICK','СПЕКТР ВХОДНОГО СИГНАЛА')

!ACCEPT*

!CALL DRAW (M,1,f,SPECTR_Ph,0,'STICK','ФАЗА ВХОДНОГО СИГНАЛА')

!ACCEPT*

! ВТОРОЕ ЗАДАНИЕ - ПОСТРОЕНИЕ МАТ. МОДЕЛИ СИСТЕМЫ.I=1,M(I)=2.0*pi*f(I)*CJ

ENDDO

!СОЗДАНИЕ МАССИВА ПЕРЕДАТОЧНЫХ ХАРАКТЕРИСТИК БЛОКОВ 9, 11, 18

DO I=1,M(I)=(CP(I)*L1+R2)/(CP(I)**2*R1*C1*L1+CP(I)*(R1*R2*C1+L1)+R1+R2)(I)=(CP(I)*L_1+R_2)/(CP(I)*L_1+R_1+R_2)(I)=(CP(I)**2*R__1*C__1*L__1+CP(I)*(R__1*R__2*C__1+L__1)+R__2)/(CP(I)**2*R__1*C__1*L__1+CP(I)*(R__1*R__2*C__1+L__1)+R__1+R__2)

END DO

!ПОЛУЧЕНИЕ АМПЛИТУДНых ХАРАКТЕРИСТИК БЛОКОВ

DO I=1,M(I)=CABS(CW9(I))(I)= CABS(CW11(I))(I)= CABS(CW18(I))

END DO

! ПОСТРОЕНИЕ ПЕРЕДАТОЧНЫХ ХАРАКТЕРИСТИК БЛОКОВ И ЧХ Т-БЛОКА

!CALL DRAW (M,1,f,W9,0,'LINE','ПЕРЕДАТОЧНАЯ ХАРАКТЕРИСТИКА БЛОКА 9')

!ACCEPT*

!CALL DRAW (M,1,f,W11,0,'LINE','ПЕРЕДАТОЧНАЯ ХАРАКТЕРИСТИКА БЛОКА 11')

!ACCEPT*

!CALL DRAW (M,1,f,W18,0,'LINE','ПЕРЕДАТОЧНАЯ ХАРАКТЕРИСТИКА БЛОКА 18')

!ACCEPT*

!CALL DRAW (K,1,FT,AChH,0,'LINE','AЧХ Т БЛОКА')

!ACCEPT*

!CALL DRAW (K,1,FT,FChH,0,'LINE','ФЧХ Т БЛОКА')

!ACCEPT*

ENDDO

! ПОСТРОЕНИЕ АПРОКСИМИРОВАНЫХ ЧХ Т-БЛОКА

!CALL DRAW (M,1,F,AChHT_EDENICY,0,'LINE','АПРОКСИМИРОВАНАЯ AЧХ Т-БЛОКА')

!ACCEPT*

!CALL DRAW (M,1,F,FChHT_RADIANY,0,'LINE','АПРОКСИМИРОВАНАЯ ФЧХ Т-БЛОКА')

!ACCEPT*

! СОЗДАНИЕ КОМПЛЕКСНОГО КОЕФИЦИЕНТА ПЕРЕДАЧИ Т-БЛОКА

DO I=1,M(I)=AChHT_EDENICY(I)*EXP(CJ*FChHT_RADIANY(I))I=1,M(I)=CABS(CWT(I))

ENDDO

!CОЗДАНИЕ КОЕФИЦИЕНТА ПЕРЕДАЧИ ВСЕЙ СИСТЕМЫ

! БЛОК 11 И Т-БЛОК ПОСЛЕДОВАТЕЛЬНОI=1,M

CW1(I)=CW11(I)*CWT(I)

ENDDO

! ПОЛУЧЕНЫЙ БЛОК 1 И БЛОК 18 С ОТРИЦАТЕЛЬНОЙ ОБР. СВ.

DO I=1,M(I)=CW1(I)/(1+CW1(I)*CW18(I))

ENDDO

! ПОЛУЧЕНЫЙ БЛОК 2 И БЛОК 9 ПОСЛЕДОВАТЕЛЬНО

DO I=1,M(I)=CW2(I)*CW9(I)

ENDDO

!МОДЕЛИРОВАНИЕ АЧХ ВСЕГО БЛОКАI=1,M

W_USTROJSTVA(I)=CABS(CW3(I))

!CALL DRAW(M,1,F,W_USTROJSTVA,0,'LINES','ОБЩАЯ АМПЛИТУДО-ЧАСТОТНАЯ ХАРАКТЕРИСТИКА ВСЕГО БЛОКА')

!ACCEPT*

! СОЗДАНИЕ СИГНАЛА НА ВЫХОДЕ СИСТЕМЫ

DO I=1,M-1(M+I+1)=CONJG(CW3(M-I+1))I=1,N(I)=CS(I)*CW3(I)FFT(CSY,N,1)I=1,N(I)=REAL(CSY(I))FFT(CSY,N,0)I=1,M(I)=REAL(CSY(I))(I)=AIMAG(CSY(I))

!CALL DRAW(N,1,t,SY,0,'LINES','СИГНАЛ НА ВЫХОДЕ')

!ACCEPT*

!CALL DRAW(M,1,F,SpRe2,0,'STIKC','РЕАЛЬНАЯ СОСТАВЛЯЮЩАЯ СИГНАЛА НА ВЫХОДЕ')

!ACCEPT*

!CALL DRAW(M,1,F,SpIm2,0,'STIKC','МНИМАЯ СОСТАВЛЯЮЩАЯ СИГНАЛА НА ВЫХОДЕ')

!ACCEPT*

!СОЗДАНИЕ МАССИВОВ АМПЛИТУДНОЙ И ФАЗОВОЙ ХАРАКТЕРИСТИК

DO I=1,M(I)=CABS(CSY(I))((ABS(REAL(CSY(I)))<1E-3).and.ABS(IMAG(CSY(I)))<1E-3) THEN(I)=0.0(I)=ATAN2(AIMAG(CSY(I)),REAL(CSY(I)))

ENDIF

!ПОСТРОЕНИЕ ЧХ ВЫХОДНОГО СИГНАЛА

!CALL DRAW(M,1,F,SpA2,0,'STIKC','АМПЛИТУДНО-ЧАСТОТНАЯ ХАРАКТЕРИСТИКА ВЫХОДНОГО СИГНАЛА')

!ACCEPT*

!CALL DRAW(M,1,F,SpPh2,0,'STIKC','ФАЗО-ЧАСТОТНАЯ ХАРАКТЕРИСТИКА ВЫХОДНОГО СИГНАЛА')

!ACCEPT*

! ТРЕТЬЕ ЗАДАНИЕ

!применение мультипликативного метода

DO I=1,N=XSTART*XKOEF=XSTART-AINT(XSTART)=XSTART(I)=-(LOG(1-R))/LYAMDA

! X(I) ПОЛУЧЕНО МЕТОДОМ ОБРАТНОЙ ФУНКЦИИ

ENDDO

!CALL DRAW(N,1,T,X,0,'LINE','ШУМ')

!ACCEPT*

! НАХОЖДЕНИЕ ДИСПЕРСИИ И СР. КВ. ОТКЛОНЕНИЯ=0.0I=1,N=DX0+(X(I)-MX)**2=DX0/REAL(N-1)=SQRT(DX)

!PRINT *,'DX=',DX,'SKO=',SIGMAX

ACCEPT*

! НАХОЖДЕНИЕ МИН. И МАКС. ЭЛЕМЕНТОВ МАССИВА=0.0=X(1)

DO I=1,N(X(I)>XMAX) XMAX=X(I)

!PRINT *,'XMIN=',XMIN, 'XMAX=',XMAX

! ПОСТРОЕНИЕ ГИСТОГРАМЫRHIST(H,NHIST,X,N,XMIN,XMAX)I=1,NHIST(I)=XMIN+REAL(I-1)*(XMAX-XMIN)/(NHIST-1)

!ACCEPT*

!PRINT *,'TT=',TT

!ACCEPT*

!CALL DRAW(NHIST,1,TT,H,0,'BOX','ГИСТОГРАММА')

!ACCEPT*

!ПОСТРОЕНИЕ ФУНКЦИИ ЭКСПОНЕНТЫI=1,NPOINT

XX(I)=XMIN+((XMAX-XMIN)/(NPOINT-1))*REAL(I-1)I=1,NPOINT(I)=LYAMDA*EXP(-LYAMDA*XX(I))

!CALL DRAW(NPOINT,1,XX,FX,0,'LINES','ЭКСПОНЕНТА')

!ACCEPT*=(XMAX-XMIN)/REAL(NPOINT)PL(PLOSHA,NHIST,INTERVAL,FX,NPOINT)I=1,NHIST(I)=PLOSHA(I)*REAL(N)=0.0I=1,NHIST=SP+PLOSHA(I)*,PLOSHA,'SP',SP*

!CALL DRAW (NHIST,1,TT,PLOSHA1,0,'LINES','ПЛОЩАДЬ')

!ACCEPT*I=1,NHIST(I,2)=H(I)(I,1)=PLOSHA1(I)(I,1)=TT(I)(I,2)=TT(I)

!CALL DRAW(NHIST, NG, TT1, DVA, (/0,0/),(/'LINES','BOX'/),'СРАВНЕНИЕ')

!ACCEPT*=0.0I=1,NHIST=XIXI+(PLOSHA1(I)-H(I))**2/PLOSHA1(I)

!PRINT *,'XIXI=', XIXI

!! 4 ЗАДАНИЕI=1,M(I)=CMPLX(X(I)-MX,0.0)FFT (CX,N,0)I=1,N(I)=CX(I)*CW3(I)

!CALL DRAW (M,1,F,Sxvyh,0,'LINES','CПЕКТР ШУМА НА ВЫХОДЕ')

!ACCEPT*I=1,M_Re(I)=REAL(Sxvyh(I))

END DO

!CALL DRAW (M,1,f,SPECTRshuma_Re,0,'STICK','РЕАЛЬНАЯ СОСТАВЛЯЮЩАЯ ВЫХОДНОГО ШУМА')

!ACCEPT*

!CALL DRAW (M,1,f,SPECTRshuma_Im,0,'STICK','МНИМАЯ СОСТАВЛЯЮЩАЯ ВЫХОДНОГО ШУМА')

!ACCEPT*I=1,M_A(I)=CABS(Sxvyh(I))

IF(ABS(REAL(Sxvyh(I)))<1.0e-3.AND.ABS(AIMAG(Sxvyh(I)))<1.0e-3) THEN

SPECTRshuma_Ph(I)=0.0

ELSE_Ph(I)=ATAN2(AIMAG(Sxvyh(I)),REAL(Sxvyh(I)))DO

!CALL DRAW (M,1,f,SPECTRshuma_A,0,'STICK','АМПЛИТУДНЫЙ СПЕКТР ВЫХОДНОГО ШУМА')

!ACCEPT*

!CALL DRAW (M,1,f,SPECTRshuma_Ph,0,'STICK','ФАЗА ВЫХОДНОГО ШУМА')

!ACCEPT*I=1,N(I)=Sxvyh(I)*CONJG(Sxvyh(I))FFT(AKF,N,1)I=1,N(I)=REAL(AKF(I))=AKF(1)I=1,N(I)=RAKF(I)/C

!CALL DRAW (N/4,1,T,RAKF1,0,'LINE','AKF')

!ACCEPT*FFT(Sxvyh,N,1)I=1,N(I)=REAL(Sxvyh(I))=SxvyhR(1)=SxvyhR(1)I=1,N(SxvyhR(I)>XMAX2) XMAX2=SxvyhR(I)(SxvyhR(I)<XMIN2) XMIN2=SxvyhR(I)*,'XMIN2=',XMIN2, 'XMAX2=',XMAX2RHIST(H1,NHIST,SxvyhR,N,XMIN2,XMAX2)I=1,NHIST(I)=XMIN2+REAL(I-1)*(XMAX2-XMIN2)/(NHIST-1)

!CALL DRAW (M,1,T,SxvyhR,0,'LINE','ШУМ НА ВЫХОДЕ СИСТЕМЫ')

!ACCEPT*

!CALL DRAW(NHIST,1,TT2,H1,0,'BOX','ГИСТОГРАММА ВЫХОДНОГО ШУМА')

!ACCEPT*

!MX0Sh=0.0I=1,NSh=MX0Sh+SxvyhR(I)=MX0Sh/REAL(N)Sh=0.0I=1,NSh=DX0Sh+(SxvyhR(I)-MXSh)**2=DX0Sh/REAL(N-1)=SQRT(DXSh)=0.0I=1,N=ALFAX0+((SxvyhR(I)-MXSh)**3)=ALFAX0/REAL(N*SIGMAXSh**3)=0.0I=1,N=EX0+((SxvyhR(I)-MXSh)**4)=EX0/(REAL(N)*(SIGMAXSh**4))-3*, 'MXSh=',MXSh, ' DXSh=',DXSh, 'ALFAX=', ALFAX, 'EX=', EX, 'SIGMAXSh=', SIGMAXSh

DO I=1,NPOINT

! ПОСТРОЕНИЕ ПЛОТНОСТИ НОРМАЛЬНОГО РАССПРЕДЕЛЕНИЯ

XNORM(I)=XMIN2+REAL(I-1)*(XMAX2-XMIN2)/REAL(NPOINT-1)=MXSh+3.0*SIGMAXSh=MXsH-3.0*SIGMAXShI=1,NPOINT(I)=(1.0/(SQRT(2.0*PI)*SIGMAXSh))*EXP(-(1.0/(2.0*SIGMAXSh**2))*(XNORM(I)-MXSh)**2)

!CALL DRAW(NPOINT,1,XNORM,PLOTN,0,'LINE','НОРМАЛЬНОЕ РАСПРЕДЕЛЕНИЕ')

!ACCEPT*

! ФОРМИРОВАНИЕ МАССИВА КОЛИЧЕСТВА ПОПАДАНИЙ В ИНТЕРВАЛ ГИСТОГРАММЫ

INTERVAL2=(XMAX2-XMIN2)/REAL(NPOINT)PL(POPADANIE,NHIST,INTERVAL2,PLOTN,NPOINT)I=1,NHIST(I)=POPADANIE(I)*REAL(N)I=1,NHIST(I,2)=H1(I)(I,1)=POPADANIE1(I)(I,1)=TT2(I)(I,2)=TT2(I)

!CALL DRAW(NHIST, NG, XNORMAL, DVA2, (/0,0/),(/'LINES','BOX'/),'СРАВНЕНИЕ')

!ACCEPT*=0.0I=1,NHIST=XIXINORM+(POPADANIE1(I)-H1(I))**2/POPADANIE1(I)*,'XIXINORM=', XIXINORM

! ПОСТРОЕНИЕ ПЛОТНОСТИ РАССПРЕДЕЛЕНИЯ РЕЛЛЕЯ

BREL=XMAX2

AREL=0.0I=1,NPOINT(I)=AREL+REAL(I-1)*(BREL-AREL)/REAL(NPOINT-1)I=1,NPOINT(I)=(XREL(I)/SIGMAXSh**2)*EXP(-XREL(I)**2/(2.0*SIGMAXSh**2))

!CALL DRAW(NPOINT,1,XREL,PLOTN2,0,'LINE','РАСПРЕДЕЛЕНИЕ РЕЛЕЯ')

!ACCEPT*

!ФОРМИРОВАНИЕ МАССИВА КОЛИЧЕСТВА ПОПАДАНИЙ В ИНТЕРВАЛ ГИСТОГРАММЫ

INTERVAL3=(BREL-AREL)/REAL(NPOINT)PL(POPADANIE2,NHIST,INTERVAL3,PLOTN2,NPOINT)I=1,NHIST(I)=POPADANIE2(I)*REAL(N)

!CALL DRAW(NHIST,1,TT2,POPADANIE2,0,'LINE','ПЛОЩАДЬ')

!ACCEPT*I=1,NHIST(I,2)=H1(I)(I,1)=POPADANIE3(I)(I,1)=TT2(I)(I,2)=TT2(I)

!CALL DRAW(NHIST, NG, XREL, DVA3, (/0,0/),(/'LINES','BOX'/),'СРАВНЕНИЕ')

!ACCEPT*=0.0I=1,NHIST=XIXIREL+(POPADANIE3(I)-H1(I))**2/POPADANIE3(I)*,'XIXIREL=', XIXIREL

ПОДПРОГРАММА ИНТЕРПОЛЯЦИИFUNCTION RINTERPOL(L, X_ZAD, f_ZAD, X)LX_ZAD(L), f_ZAD(L), XCH(L),ZN(L),WTI=1,L(I)=1.0(I)=1.0J=1,L(I/=J) CH(I)=CH(I)*(X-X_ZAD(J))(I/=J) THEN(I)=ZN(I)*(X_ZAD(I)-X_ZAD(J))=0.0I=1,L=WT+f_ZAD(I)*CH(I)/ZN(I)=WTFUNCTION RINTERPOL

ПОДПРОГРАММА ПОДСЧЕТА ПЛОЩАДИ

SUBROUTINE PL(S,N1,H,Y,N)

INTEGER N,N1Y(N)SUMMA,H,S(N1)J=2,N-1,2(J)=4.0*Y(J)I=1,N,2(I)=2.0*Y(I)J=1,N1=0.0I=(J-1)*INT(N/N1)+1,J*INT(N/N1)+1((J-1)*INT(N/N1)+1)=Y((J-1)*INT(N/N1)+1)/2.0(J*INT(N/N1)+1)=Y(J*INT(N/N1)+1)/2.0=SUMMA+Y(I)(J)=SUMMA*H/3.0

ENDDO