# Функция для парсинга одной строки
function parse_line(line)
    # Заменяем запятые на точки для корректного парсинга чисел
    line = replace(line, "," => ".")
    # Разделяем по пробелам/табам и фильтруем пустые элементы
    elements = filter(!isempty, split(line, r"\s+"))
    # Преобразуем элементы в нужные типы
    co = parse(Float64, elements[1])
    ch = parse(Float64, elements[2])
    co2 = parse(Float64, elements[3])
    o2 = parse(Float64, elements[4])
    lamb = parse(Float64, elements[5])
    n = parse(Int, elements[6])
    nox = parse(Int, elements[7])
    tmas = parse(Float64, elements[8])
    time = Time(elements[9])
    return (CO=co, CH=ch, CO2=co2, O2=o2, lamb=lamb, n=n, NOx=nox, Tmas=tmas, Time=time)
end

using DataFrames
# Чтение данных
data_lines = readlines("$(@__DIR__)/data.txt")[2:end]  # Пропускаем заголовок
parsed_data = [parse_line(line) for line in data_lines] # Парсим все строки
df = DataFrame(parsed_data) # Преобразуем в DataFrame
println(first(df, 5)) # Выводим первые 5 строк для проверки
println()
using Dates, FFTW, Statistics
# Добавим колонку с временем в секундах
start_time = df.Time[1]
df[!, :Seconds] = [Dates.value(t - start_time) / 1e9 for t in df.Time]  # переводим наносекунды в секунды
println(first(df, 5)) # Выводим первые 5 строк для проверки

i = 1 # Номер столбца для анализа
y = df[:,i]
column_names = names(df) # Получаем имена столбцов
println("Имя столбца для анализа: $(column_names[i])")
t = df.Seconds.-(15*60*60+56*60); t[1] = 0;
println("t: $(t[1:5])") # Выводим первые 5 строк для проверки

window_size = 15  # размер окна для сглаживания
trend = [mean(y[max(1, i-window_size):min(end, i+window_size)]) for i in 1:length(y)]

# График исходных данных и тренда
plot(t, y, label="Исходные данные", xlabel="Время (с)", ylabel="CO", legend=:topleft)
plot!(t, trend, label="Тренд", linewidth=2)

using Polynomials
trend_coef = fit(t, y, 3)  # полином третьей степени
trend = trend_coef.(t)

# График исходных данных и тренда
plot(t, y, label="Исходные данные", xlabel="Время (с)", ylabel="CO", legend=:topleft)
plot!(t, trend, label="Тренд", linewidth=2)

# Удаляем тренд для анализа периодических составляющих
detrended = y .- trend

# Вычисляем БПФ
n = length(detrended)
freq = fftfreq(n, 1/mean(diff(t)))  # частоты в Гц
fft_vals = fft(detrended)

# Амплитуды и фазы
amplitude = abs.(fft_vals) ./ n * 2
phase = angle.(fft_vals)

# Оставляем только положительные частоты (симметрия БПФ)
pos_freq = freq[1:div(n,2)+1]
pos_ampl = amplitude[1:div(n,2)+1]

# График амплитудного спектра
plot(pos_freq, pos_ampl, xlabel="Частота (Гц)", ylabel="Амплитуда", 
     title="Амплитудный спектр", legend=false)

# Находим 3 наиболее значимые частоты
top_n = 3
sorted_idx = sortperm(pos_ampl, rev=true)
top_freq = pos_freq[sorted_idx[1:top_n]]
top_ampl = pos_ampl[sorted_idx[1:top_n]]
top_phase = phase[sorted_idx[1:top_n]]

# Восстанавливаем гармоники
harmonics = zeros(n)
for (f, a, ϕ) in zip(top_freq, top_ampl, top_phase)
    harmonics .+= a .* cos.(2π * f * t .+ ϕ)
end

# Итоговый график
plot(t, detrended, label="Детрендированные данные", xlabel="Время (с)", ylabel="CO")
plot!(t, harmonics, label="Основные гармоники", linewidth=2)

plot(t, y, label="Исходные данные", xlabel="Время (с)", ylabel="CO")
plot!(t, trend, label="Тренд", linewidth=2)
plot!(t, trend .+ harmonics, label="Тренд + гармоники", linewidth=2, linestyle=:dash)