%% 生成基本表格信息
%% 加载路径
fun_path = genpath('functions');
addpath(fun_path)
%% 基本信息
region = 'US'; % 地区
Root_Path1 = ['E:\2024.10 甲烷比赛数据\', region]; % 原始数据
Root_Path2 = 'E:\2024.10 甲烷比赛结果'; % 核查结果
% Zone = 'MST'; % 地区时 亚利桑那
Zone = 'CST'; % 德克萨斯
%% 获取文件夹信息
folders = dir(Root_Path1);
folders = folders([folders.isdir]);
folders = folders(~ismember({folders.name}, {'.', '..'}));
len = size(folders, 1);
% 初始化表格
tableData = cell(9999, 17);
%% 获取已知点源坐标
tablename = 'data/点源位置/点源坐标.xlsx';
T_coord = readtable(tablename, 'VariableNamingRule', 'preserve', 'TreatAsEmpty', '');
len_coord = size(T_coord, 1); % 数据库中 点源的个数
%% 循环处理文件夹
k = 0;
for i = 1 : len
    %% 数据名称
    filename = folders(i).name;
    origin_path = fullfile(Root_Path1, filename);
    flist = ReadFiles(origin_path); % 读取数据
    xmlStruct = ReadXmlData(flist.xml); % 读取xml文件
    %% 获取卫星名称
    satellite = xmlStruct.SatelliteID;
    %% 获取8位数日期
    underIndices = strfind(filename, '_');
    date = filename(underIndices(4) + 1 : underIndices(4) + 8);
    %% 创建景文件夹名称
    dataname = filename(underIndices(end) + 1 : end);
    foldername = [region, '_', date, '_', dataname];
    %%
    [ToS, UTC, LST] = processTime(xmlStruct, Zone);
    %% 判断该景数据是否覆盖已知源
    for j = 1:len_coord
        facility_name = char(T_coord{j, 'name'});
        type = char(T_coord{j, 'type'});
        lon = double(T_coord{j, 'lon'});
        lat = double(T_coord{j, 'lat'});
        isInside = isPointInTiffRegion(lon, lat, xmlStruct);
        if(isInside)
            subfoldername = facility_name;
            %% 保存到表格
            k = k + 1;
            tableData{k, 1} = k; % 序号
            tableData{k, 2} = filename; % 文件名
            tableData{k, 3} = foldername;
            tableData{k, 4} = facility_name; % 设施名称
            tableData{k, 5} = type;
            tableData{k, 6} = lat;
            tableData{k, 7} = lon;
            tableData{k, 8} = satellite;
            tableData{k, 9} = date; % 8位日期
            tableData{k, 10} = ToS;
            tableData{k, 11} = UTC;
            tableData{k, 12} = LST;
        end
    end
end
tableData = tableData(1:k, :);
%% 创建表格
folderTable = cell2table(tableData, 'VariableNames', ...
    {'id', 'filename', 'foldername', 'facility', 'type', 'Latitude', 'Longitude', ...
    'satellite', 'date', 'Time of Satellite', 'UTC', 'Local Standard Time', ...
    '反演结果', '不确定度', '风速', '排放', '备注'});
% 显示表格
disp(folderTable);
%% 导出表格
outputname = ['data/', '2024.10 甲烷比赛', region, '.xlsx'];
writetable(folderTable, outputname);
%% 创建文件夹
folderNames = tableData(:, 3);
subfolderNames = tableData(:, 4);
% 遍历文件夹名称并创建子文件夹
for i = 1:length(folderNames)
    FolderPath = fullfile(Root_Path2, folderNames{i});
    mymkdir(FolderPath);
    SubFolderPath = fullfile(Root_Path2, folderNames{i}, subfolderNames{i});
    mymkdir(SubFolderPath); % 创建路径
end
disp('文件夹创建完成。');

%% 判断是否包含在内
function isInside = isPointInTiffRegion(lon, lat, xmlStruct)
% Extract the boundary coordinates from xmlStruct
topLeftLat = str2double(xmlStruct.TopLeftLatitude);
topLeftLon = str2double(xmlStruct.TopLeftLongitude);
topRightLat = str2double(xmlStruct.TopRightLatitude);
topRightLon = str2double(xmlStruct.TopRightLongitude);
bottomLeftLat = str2double(xmlStruct.BottomLeftLatitude);
bottomLeftLon = str2double(xmlStruct.BottomLeftLongitude);
bottomRightLat = str2double(xmlStruct.BottomRightLatitude);
bottomRightLon = str2double(xmlStruct.BottomRightLongitude);

% Define the coordinates of the corners of the TIFF region
lons = [topLeftLon, topRightLon, bottomRightLon, bottomLeftLon, topLeftLon];
lats = [topLeftLat, topRightLat, bottomRightLat, bottomLeftLat, topLeftLat];

% Check if the point (lon, lat) is within the polygon defined by the corners
isInside = inpolygon(lon, lat, lons, lats);

% 绘制TIFF区域和点
% figure;
% hold on;
% plot(lons, lats, '-o', 'LineWidth', 2, 'DisplayName', 'TIFF Region'); % 绘制多边形
% if isInside
%     plot(lon, lat, 'g*', 'MarkerSize', 10, 'DisplayName', 'Point (Inside)'); % 如果点在多边形内，使用绿色星形标记
% else
%     plot(lon, lat, 'r*', 'MarkerSize', 10, 'DisplayName', 'Point (Outside)'); % 如果点在多边形外，使用红色星形标记
% end
% xlabel('Longitude');
% ylabel('Latitude');
% title('TIFF Region and Point');
% legend;
% grid on;
% hold off;
end

%% 时间
function [ToS, UTC, LST] = processTime(xmlStruct, Zone)
%% 卫星时间
ToS = datetime(xmlStruct.StartTime, 'Format', 'yyyy-MM-dd HH:mm:ss'); % 卫星时间
%% UTC时间
switch xmlStruct.SatelliteID
    case 'GF5B'
        UTC = ToS - duration(8, 0, 0);
    case 'ZY1E'
        UTC = ToS;
    case 'ZY1F'
        UTC = ToS;
    case 'GF5'
        UTC = ToS - duration(8, 0, 0);
    case 'GF5A'
        UTC = ToS;
end
%% 当地时间
switch Zone
    case 'EST'
        LST = UTC - duration(5, 0, 0);
    case 'CST'
        LST = UTC - duration(6, 0, 0);
    case 'MST'
        LST = UTC - duration(7, 0, 0);
    case 'PST'
        LST = UTC - duration(8, 0, 0);
    case 'Beijing'
        LST = UTC + duration(8, 0, 0); % 北京时间
    case 'Warsaw'
        LST = UTC + duration(2, 0, 0); % 华沙时间
    otherwise
        disp('ERR Zone name!')
end
end