Merge branch '1-matlab-scripts' into 'main'

add Matlab_Scripts

Closes #1

See merge request !4

Matlab_Scripts/ReadMeMatlab.dat

+1. Необходимо установить себе на компьютер Matlab (если его нет).
+
+2. В Matlabe уже находится огромное количество встроенных функций. Однако для некоторых целей их может не хватать. Например, для
+отрисовки данных на карте необходима внешняя функция m_map (есть и встроенные функции, но с этой работать удобнее).
+Её можно самостоятельно встроить в Matlab. Для этого необходимо разархивировать папку m_map и переместить её в папку C/.../MATLAB/toolbox. 
+Затем отредактировать файл C/.../MATLAB/toolbox/local/pathdef.m, вставив туда по образцу строчку с путем к папке m_map. Скопируйте файл 
+pathdef.m к себе на рабочий стол. Вставьте необходимую строку (путь к папке m_map). И переместите файл pathdef.m обратно в папку local 
+(с перезаписью старого файла). Далее, перезапустив Matlab, ввести команду rehash toolboxcache. На всякий случай можно перезапустить Matlab.
+
+3. Код программы обычно пишется в так называемом Script (Скрипте). Перед комментариями ставится знак процента (%). Таким образом Matlab 
+понимает, что это не команда. Можно разделять Скрипт на блоки, ставя два знака процента (%%) подряд.
+Чтобы запустить весь Скрипт, нужно нажать на кнопку Run.
+Если вы хотите запустить только какую-то часть (блок), нужно поставить курсор в любую часть этого блока и нажать Run Section (или 
+комбинацию клавиш Cntrl+Enter).
+Также можно запускать (отправлять в командную строку) любую часть Скрипта (не привязываясь к конкретным блокам). Для этого нужно
+выделить кусок программы, который вы хотите запустить, и нажать f9. 
+
+4. Для ознамления с возможностями Matlab'a даны два тестовых Скрипта  - Test1.m и Test2.m. Обязательно прописать в них путь к папке, куда
+вы положите тестовые данные реанализа!
+Успехов!
+

Matlab_Scripts/Test1.m

+% Очищаем командную строку (Command Window)
+clc
+% Очищаем Workspace (удаляем все переменные)
+clear all
+
+%% Открываем данные
+% Прописываем путь к файлу
+% Знак точка с запятой ";" в конце строки блокирует вывод полученного в ходе 
+% выполнения команды результата в командной строке
+filepath = "C:\";
+% Прописываем название файла
+filename = "surface_t2m_precip.nc";
+% Просматриваем содержимое файла NetCDF командой ncdisp (переменные, их размерность и др.)
+% Информация отображается в командной строке 
+ncdisp (filepath + filename);
+
+% Выбираем необходимые нам переменные
+% Загружаем в Workspace температуру на 2 метрах командой ncread
+% В скобках указываем название переменной - то, как она называется в файле
+% В Workspace должна появиться переменная t2m_data в виде массива
+% размерностью 1440*721*3 (долготы*широты*время)
+t2m_data = ncread(filepath + filename, 't2m');
+% Аналогичным образом загружаем векторы с долготами, широтами и временем
+% Если нажать на появившиеся перемнные, в окошке Variables можно
+% посмотреть, что содержится внутри них
+lonn = ncread(filepath + filename, 'longitude');
+latt = ncread(filepath + filename, 'latitude');
+time = ncread(filepath + filename, 'time');
+
+%% Преобразовываем данные
+% Данные по температуре даны в Кельвинах -> переводим их в градусы Цельсия
+t2m_data_Celsius = t2m_data - 273.15;
+
+% Данные о времени даны в количестве часов с 00:00 1 января 1900 г. -> переводим их в понятные нам даты
+% В Matlab'е есть команда datenum для работы с данными формата даты
+% Она переводит заданную дату в количество дней с 00:00 1 января 0000 г.
+% создадим переменную t0 - количество дней с 00:00 1 января 0000 г. по 00:00 1 января 1900 г.
+t0 = datenum('1900-01-01 00:00:00', 'yyyy-mm-dd HH:MM:SS');
+% переведем даты из нашего файла из формата "часы с 00:00 1 января 1900 г."
+% в "дни с с 00:00 1 января 1900 г."
+% time_days_since1900year - количество дней с 00:00 1 января 1900 г.по
+% даты, указанные в файле. Точка перед знаком деления обозначает, что
+% данная арифметическая операция применяется отдельно к каждому члену
+% массива
+time_days_since1900year = time ./ 24;
+% сложим две переменные между собой, чтобы получить time_days_since0year - количество дней с 00:00 1 января 0000 г. 
+% по даты, указанные в файле
+time_days_since0year = time_days_since1900year + t0;
+% конвертируем с помощью команды datetime из формата "дни с с 00:00 1 января 0000 г." в нормальный вид, 
+% указав при этом, что стартуем именно с 0000 года флагом 'datenum'
+% если нажать на переменную time_real, можно увидеть, за какие даты даны
+% данные
+time_I_understand = datetime(time_days_since0year, 'ConvertFrom', 'datenum');
+
+%% Отрисовка данных на карте
+% Переворачиваем командой flipud вектор широт, чтобы они начинались с -90
+latt1 = flipud(latt);
+% Строим координатную сетку с размерностью долгот и широт
+[lon,lat] = meshgrid(lonn, latt1);
+
+% Открываем окно для нашего рисунка во весь экран
+figure('units','normalized','outerposition',[0 0 1 1])
+% Указываем проекцию и координаты
+m_proj ('Equidistant Cylindrical','lon',[0 360], 'lat', [-90 90]);
+% Выбираем переменную за первый шаг по времени 
+t2m_data_one_time = t2m_data_Celsius(:,:,1);
+% Транспонируем (знак апострофа "'") -> меняем местами
+% размерности), так как сетка имеет размерность широты*долготы
+t2m_data_Celsius_transp = t2m_data_one_time';
+% Отрисовываем переменную, заодно переворачивая её командой flipud, чтобы
+% данные начинались с 90° ю.ш.
+m_contourf(lon, lat, flipud(t2m_data_Celsius_transp),'LineColor','none', 'LevelStep', 5);
+% Отрисовываем линию побережья
+m_coast('Color', 'k','LineWidth', 1)
+% Наносим линии сетки
+m_grid
+% Обозначаем цветовую шкалу
+colormap(turbo)
+% Выводим цветовую шкалу на рисунок
+c = colorbar;
+% Подписываем шкалу
+c.Label.String = 'Temperature, °C';
+% Поворачиваем её на -90 градусов
+c.Label.Rotation = -90;
+% Сдвигаем её, чтобы она не заезжала на шкалу после поворота
+c.Label.Position = [3.5 -12.5 0];
+% Определяем размер шрифта подписи
+c.Label.FontSize = 12;
+
+% Добавляем название рисунка
+% Команада datestr переводит значение переменной из формата даты в формат
+% строки
+% По умолчанию название будет написано жирным шрифтом -> команда 'FontWeight','Normal'
+% делает подпись обычной
+title("Air temperature at 2 m - " + datestr(time_I_understand(1), 'mmmm, yyyy'),'FontWeight','Normal','FontSize',16)
+    
+%% Отрисовка данных на карте циклом для всех точек по времени и сохранение рисунков
+% Цикл от 1 до длины вектора с точками по времени, то есть от 1 до 3 с
+% шагом 1
+for i = 1:length(time_I_understand)
+    % Открываем окно для нашего рисунка во весь экран
+    figure('units','normalized','outerposition',[0 0 1 1])
+    % Указываем проекцию и координаты
+    m_proj ('Equidistant Cylindrical','lon',[0 360], 'lat', [-90 90]);
+    % Выбираем переменную за конкретный шаг по времени 
+    t2m_data_one_time = t2m_data_Celsius(:,:,i);
+    % Транспонируем (знак апострофа "'") -> меняем местами
+    % размерности), так как сетка имеет размерность широты*долготы
+    t2m_data_Celsius_transp = t2m_data_one_time';
+    % Отрисовываем переменную, заодно переворачивая её командой flipud, чтобы
+    % данные начинались с 90° ю.ш.
+    m_contourf(lon, lat, flipud(t2m_data_Celsius_transp),'LineColor','none', 'LevelStep', 5);
+    % Отрисовываем линию побережья
+    m_coast('Color', 'k','LineWidth', 1)
+    % Наносим линии сетки
+    m_grid
+    % Обозначаем цветовую шкалу
+    colormap(turbo)
+    % Выводим цветовую шкалу на рисунок
+    c = colorbar;
+    % Подписываем шкалу
+    c.Label.String = 'Temperature, °C';
+    % Поворачиваем её на -90 градусов
+    c.Label.Rotation = -90;
+    % Сдвигаем её, чтобы она не заезжала на шкалу после поворота
+    c.Label.Position = [3.5 -12.5 0];
+    % Определяем размер шрифта подписи
+    c.Label.FontSize = 12;
+    % Ограничиываем размах значений шкалы
+    clim([-70 45])
+    % Добавляем название рисунка
+    % Команада datestr переводит значение переменной из формата даты в формат
+    % строки
+    % По умолчанию название будет написано жирным шрифтом -> команда 'FontWeight','Normal'
+    % делает подпись обычной
+    title("Air temperature at 2 m - " + datestr(time_I_understand(i), 'mmmm, yyyy'),'FontWeight','Normal','FontSize',16)
+    
+    % Сохраняем рисунок в формате png с разрешением 300 dpi
+    print("T2m_temperature_" + datestr(time_I_understand(i), 'mmmm_yyyy'), '-dpng', '-r300');
+
+end
+
+
+
+
+

Matlab_Scripts/Test2.m

+% Очищаем командную строку (Command Window)
+clc
+% Очищаем Workspace (удаляем все переменные)
+clear all
+
+%% Открываем данные
+% Прописываем путь к файлу
+% Знак точка с запятой ";" в конце строки блокирует вывод полученного в ходе 
+% выполнения команды результата в командной строке
+filepath = "C:\";
+% Прописываем название файла
+filename = "p-levels_t_u_v_geopotential.nc";
+% Просматриваем содержимое файла NetCDF командой ncdisp (переменные, их размерность и др.)
+% Информация отображается в командной строке 
+ncdisp (filepath + filename);
+
+% Выбираем необходимые нам переменные
+% Загружаем в Workspace u- и v- компоненты скорости ветра  командой ncread
+% В скобках указываем название переменной - то, как она называется в файле
+% В Workspace должны появиться переменные u_data и v_data в виде 4-х мерного массива
+u_data = ncread(filepath + filename, 'u');
+v_data = ncread(filepath + filename, 'v');
+
+% Можно посмотреть размер массива командой size и вывести в переменную u_dimensions
+% Размер массива u_data: 11440*721*3*3 (долготы*широты*барические уровни*время)
+u_dimensions = size(u_data);
+
+% Аналогичным образом загружаем векторы с долготами, широтами и временем
+% А также барические уровни (700, 850 и 1000 гПа)
+lonn = ncread(filepath + filename, 'longitude');
+latt = ncread(filepath + filename, 'latitude');
+time = ncread(filepath + filename, 'time');
+p_level = ncread(filepath + filename, 'level');
+
+%% Преобразовываем данные
+% Вычисляем полный вектор скорости через u- и v- компоненты
+wind = sqrt(u_data .^ 2 + v_data .^ 2);
+% Находим среднее значение скорости ветра за три летних месяца командой
+% mean. В скобках указывается, к какой размерности массива применяется
+% операция. В нашем случае для времени - это 4. 
+wind_ave_time = mean(wind, 4);
+
+% Данные о времени даны в количестве часов с 00:00 1 января 1900 г. -> переводим их в понятные нам даты
+% В Matlab'е есть команда datenum для работы с данными формата даты
+% Она переводит заданную дату в количество дней с 00:00 1 января 0000 г.
+% создадим переменную t0 - количество дней с 00:00 1 января 0000 г. по 00:00 1 января 1900 г.
+t0 = datenum('1900-01-01 00:00:00', 'yyyy-mm-dd HH:MM:SS');
+% переведем даты из нашего файла из формата "часы с 00:00 1 января 1900 г."
+% в "дни с с 00:00 1 января 1900 г."
+% time_days_since1900year - количество дней с 00:00 1 января 1900 г.по
+% даты, указанные в файле. Точка перед знаком деления обозначает, что
+% данная арифметическая операция применяется отдельно к каждому члену
+% массива
+time_days_since1900year = time ./ 24;
+% сложим две переменные между собой, чтобы получить time_days_since0year - количество дней с 00:00 1 января 0000 г. 
+% по даты, указанные в файле
+time_days_since0year = time_days_since1900year + t0;
+% конвертируем с помощью команды datetime из формата "дни с с 00:00 1 января 0000 г." в нормальный вид, 
+% указав при этом, что стартуем именно с 0000 года флагом 'datenum'
+% если нажать на переменную time_real, можно увидеть, за какие даты даны
+% данные
+time_I_understand = datetime(time_days_since0year, 'ConvertFrom', 'datenum');
+
+%% Отрисовка данных на карте
+% Переворачиваем командой flipud вектор широт, чтобы они начинались с -90
+latt1 = flipud(latt);
+% Строим координатную сетку с размерностью долгот и широт
+[lon,lat] = meshgrid(lonn, latt1);
+
+% Открываем окно для нашего рисунка во весь экран
+figure('units','normalized','outerposition',[0 0 1 1])
+% Указываем проекцию и координаты
+m_proj ('Equidistant Cylindrical','lon',[0 360], 'lat', [-90 90]);
+% Выбираем переменную на первом барическом уровне (700 гПа) 
+wind_one_lev = wind_ave_time(:,:,1);
+% Транспонируем (знак апострофа "'") -> меняем местами
+% размерности), так как сетка имеет размерность широты*долготы
+wind_transp = wind_one_lev';
+% Отрисовываем переменную, заодно переворачивая её командой flipud, чтобы
+% данные начинались с 90° ю.ш.
+m_contourf(lon, lat, flipud(wind_transp),'LineColor','none', 'LevelStep', 1);
+% Отрисовываем линию побережья
+m_coast('Color', 'k','LineWidth', 1)
+% Наносим линии сетки
+m_grid
+% Обозначаем цветовую шкалу
+colormap(turbo)
+% Выводим цветовую шкалу на рисунок
+c = colorbar;
+% Подписываем шкалу
+c.Label.String = 'Wind, m/s';
+% Поворачиваем её на -90 градусов
+c.Label.Rotation = -90;
+% Сдвигаем её, чтобы она не заезжала на шкалу после поворота
+c.Label.Position = [3.5 -12.5 0];
+% Определяем размер шрифта подписи
+c.Label.FontSize = 12;
+% Ограничиываем размах значений шкалы
+clim([0 25])
+% Добавляем название рисунка
+% Команада datestr переводит значение переменной из формата даты в формат
+% строки
+% По умолчанию название будет написано жирным шрифтом -> команда 'FontWeight','Normal'
+% делает подпись обычной
+title("Mean Wind Speed" + " at " + int2str(p_level(1)) + "hPa; " + datestr(time_I_understand(1), 'mmmm') + " - " + datestr(time_I_understand(3), 'mmmm, yyyy') ,'FontWeight','Normal','FontSize',16)
+
+%% Отрисовка данных на карте циклом для всех высот и сохранение рисунков
+% Цикл от 1 до длины вектора с точками по времени, то есть от 1 до 3 с
+% шагом 1
+for i = 1:length(p_level)
+    % Открываем окно для нашего рисунка во весь экран
+    figure('units','normalized','outerposition',[0 0 1 1])
+    % Указываем проекцию и координаты
+    m_proj ('Equidistant Cylindrical','lon',[0 360], 'lat', [-90 90]);
+    % Выбираем переменную на первом барическом уровне (700 гПа) 
+    wind_one_lev = wind_ave_time(:,:,i);
+    % Транспонируем (знак апострофа "'") -> меняем местами
+    % размерности), так как сетка имеет размерность широты*долготы
+    wind_transp = wind_one_lev';
+    % Отрисовываем переменную, заодно переворачивая её командой flipud, чтобы
+    % данные начинались с 90° ю.ш.
+    m_contourf(lon, lat, flipud(wind_transp),'LineColor','none', 'LevelStep', 1);
+    % Отрисовываем линию побережья
+    m_coast('Color', 'k','LineWidth', 1)
+    % Наносим линии сетки
+    m_grid
+    % Обозначаем цветовую шкалу
+    colormap(turbo)
+    % Выводим цветовую шкалу на рисунок
+    c = colorbar;
+    % Подписываем шкалу
+    c.Label.String = 'Wind, m/s';
+    % Поворачиваем её на -90 градусов
+    c.Label.Rotation = -90;
+    % Сдвигаем её, чтобы она не заезжала на шкалу после поворота
+    c.Label.Position = [3.5 -12.5 0];
+    % Определяем размер шрифта подписи
+    c.Label.FontSize = 12;
+    % Ограничиываем размах значений шкалы
+    clim([0 25])
+    % Добавляем название рисунка
+    % Команада datestr переводит значение переменной из формата даты в формат
+    % строки
+    % По умолчанию название будет написано жирным шрифтом -> команда 'FontWeight','Normal'
+    % делает подпись обычной
+    title("Mean Wind Speed" + " at " + int2str(p_level(i)) + "hPa; " + datestr(time_I_understand(1), 'mmmm') + " - " + datestr(time_I_understand(3), 'mmmm, yyyy') ,'FontWeight','Normal','FontSize',16)
+    
+    % Сохраняем рисунок в формате png с разрешением 300 dpi
+    print("Mean_Wind_Speed_at" + int2str(p_level(i)) + "hPa_" + datestr(time_I_understand(1), 'mmmm') + "-" + datestr(time_I_understand(3), 'mmmm_yyyy'), '-dpng', '-r300');
+end
+
+
+%% Строим вертикальный разрез скорости ветра, осредненной по всем широтам, долготам и за весь временной период
+% По времени мы уже осредняли, необходимо повторить операцию осреднения применительно к долготам и широтам
+% Осредняем по всем долготам. Теперь у нас не 1440 значений по долготе, а 1. Размер массива 1*721*3.
+% Команда squeeze убирает размерность, которая равна 1 => массив становится размером 721*3 
+wind_ave_lon = squeeze(mean(wind_ave_time, 1));
+% Аналогично осредняем по всем широтам
+wind_ave_lat = squeeze(mean(wind_ave_lon, 1));
+
+%% Отрисовывем полученные значения на графике командой plot
+% Указываем, что будет по оси X - скорость ветра, и что по оси Y -
+% барические уровни
+plot(wind_ave_lat, p_level, 'Color', [190/256 0 0], 'LineWidth', 2)
+% Переворачиваем ось Y, чтобы значение барических уровней шли по убыванию
+% (от 1000 к 700 гПа)
+set(gca, 'YDir','reverse')
+% Подписываем оси
+xlabel('Wind speed, m/s')
+ylabel('Pressure level, hPa')
+% Накладываем сетку
+grid on
+% Ограничиваем размах значений по оси X
+xlim([0 7])
+% Добавляем название к графику
+title("Mean Wind Speed for " + datestr(time_I_understand(1), 'mmmm') + " - " + datestr(time_I_understand(3), 'mmmm, yyyy') ,'FontWeight','Normal','FontSize',12)
+
+% Сохраняем график
+print("Vertical_plot_Mean_Wind_Speed_for_" + datestr(time_I_understand(1), 'mmmm') + "-" + datestr(time_I_understand(3), 'mmmm_yyyy'), '-dpng', '-r300');
+
+
+
+
+
+
+