%% 讀取txt文件中的站點(diǎn)顷窒，溫度畴博，降水,濕度等數(shù)據(jù)
% By Yang，2019/5/12
% 所有站點(diǎn)睛驳，所有數(shù)據(jù)都在一個(gè)txt中
% 先利用fopen打開(kāi)文件蒜田，賦予文件代號(hào)fid
fid = fopen('I:\climate\raw_climate\氣溫降水濕度1961-2018\2018new.csv');
% 利用textscan讀取txt里的數(shù)據(jù)
% %q為讀取雙引號(hào)括起來(lái)的數(shù)據(jù),共九列
% 'Delimiter',',' 以,為分隔符
% 'headerlines', 1 不讀第一行
data = textscan(fid,'%q%q%q%q%q%q%q%q','Delimiter',',','headerlines', 1);
fclose(fid);

temp1 = data{1,1}; %站點(diǎn)號(hào)V01000
tempp1 = str2double(temp1); %轉(zhuǎn)換字符串為數(shù)值
temp2 = data{1,2}; %年V04001
tempp2 = str2double(temp2);
temp3 = data{1,3}; %月V04002
tempp3 = str2double(temp3);
temp4 = data{1,4}; %日V04003
tempp4 = str2double(temp4);
temp5 = data{1,6}; %降水?dāng)?shù)據(jù)V13201稿械，溫度數(shù)據(jù)V12001,平均相對(duì)濕度數(shù)據(jù)V13003...
tempp5 = str2double(temp5);

begin_year=2018; %讀取數(shù)據(jù)的起始年份和終止年份
end_year=2018;
%組合矩陣
temm = [tempp1,tempp2,tempp3,tempp4,tempp5];

%% 按站點(diǎn)分成不同的矩陣
station = unique(temm(:,1) ); %提取不重復(fù)站點(diǎn)號(hào)
%不同站點(diǎn)的數(shù)據(jù)分為不同元胞，排為一行
for i = 1:length(station)
index = find(temm(:,1)==station(i));%利用find找到對(duì)應(yīng)站點(diǎn)數(shù)據(jù)的行號(hào)
data_d{1,i} = temm(index,:);
end

%% 記錄數(shù)據(jù)缺少太多的站點(diǎn)位置冲粤，并刪掉
%創(chuàng)建一個(gè)完整的日歷持續(xù)時(shí)間數(shù)組美莫；caldays表示間隔天數(shù)
t = datetime(begin_year,01,01):caldays(1):datetime(end_year,12,31);
[y,m,d] = ymd(t); %返回t中日期時(shí)間值的年、月和天數(shù)
d_max = length(d); %總天數(shù)(數(shù)據(jù)長(zhǎng)度)

B = zeros(1,length(data_d)); %創(chuàng)建相同數(shù)據(jù)長(zhǎng)度的0矩陣
for j = 1:length(data_d)
number = length(data_d{1,j}); %每個(gè)站點(diǎn)元胞內(nèi)的數(shù)據(jù)長(zhǎng)度
if number>round(d_max*0.9) %準(zhǔn)備把數(shù)據(jù)中缺少5%以上的刪去
B(1,j)=1; %將沒(méi)有缺少的位置標(biāo)記為1
end
end
[p,q]=find(B==0); %[p,q]為空元胞所在行列號(hào)色解，即缺少5%以上的數(shù)據(jù)
data_d(q)=[]; %去掉空元胞(即刪去),這里考慮的是原始數(shù)據(jù)是1行n列的情況，如果是m行1列就改成A(p)=[];
station(q)=[]; %對(duì)應(yīng)的站點(diǎn)號(hào)也刪掉

%% 找到每個(gè)站點(diǎn)數(shù)據(jù)中的缺測(cè)值餐茵，并填補(bǔ)為NAN
%因?yàn)閿?shù)據(jù)源有問(wèn)題科阎，缺測(cè)的值沒(méi)有填充，被直接刪除了忿族，所以要先找到缺測(cè)的天數(shù)（位置）
perfect = [y',y',m',d']; %組合為一個(gè)矩陣,第一列隨便填充锣笨，后三列年，月道批，日
for k = 1:length(data_d) %刪掉缺測(cè)站點(diǎn)所剩的站點(diǎn)數(shù)
data = data_d{1,k};
a = data_d{1,k}(1,1); %a為對(duì)應(yīng)站點(diǎn)的站點(diǎn)號(hào)
perfect(:,1) = a; %第一列替換為對(duì)應(yīng)站點(diǎn)的站點(diǎn)號(hào)错英，年月日不變
perfect(:,5) = NaN; %第五列填充為NaN
for i =1:length(data)
year = data(i,2);
month = data(i,3);
day = data(i,4);
%利用find找到年月日都正確的站點(diǎn)數(shù)據(jù)的行號(hào)
ind = find(perfect(:,2)==year&perfect(:,3)==month&perfect(:,4)==day);
perfect(ind,5) = data(i,5); %缺測(cè)的數(shù)據(jù)會(huì)沒(méi)有填充即為NaN
end
data_d2{1,k} = perfect;
end

%% 處理特征值
%降雨32700為微量，將其設(shè)為0隆豹；32744椭岩，32766等為缺測(cè)，將其轉(zhuǎn)為NAN值
for i = 1:length(data_d2)
data_d2{1,i}(data_d2{1,i} == 999999) = NaN; %溫度等數(shù)據(jù)則不需考慮該步
end

%% 將數(shù)據(jù)中NAN值較多的站點(diǎn)刪去
%統(tǒng)計(jì)每個(gè)站點(diǎn)元胞中的NaN值數(shù)量
for k = 1:length(data_d2)
temp = data_d2{1,k};
index = isnan(temp(:,5)); %利用isnan來(lái)確定是否為NaN，是為1判哥，否為0
ID = find(index == 1); %利用find找到NaN的行號(hào)
ID_temp{1,k} = temp(ID,:); %ID_temp為每個(gè)NAN值的數(shù)量
end

% 將數(shù)據(jù)中NaN值太多的站點(diǎn)刪掉
[m,n] = size(ID_temp);
B = zeros(1,n);
for j=1:n
number = length(ID_temp{1,j});
if number<round(d_max*0.05) %準(zhǔn)備把數(shù)據(jù)中缺少5%以上的刪去
B(1,j)=1; %將超過(guò)5%以上NaN的數(shù)據(jù)位置標(biāo)記為1
end
end
[p,q]=find(B==0); %[p,q]為空元胞所在行列號(hào)献雅，即缺少5%以上的數(shù)據(jù)
data_d2(q)=[]; %去掉空元胞（即刪去）,這里考慮的是原始數(shù)據(jù)是1行n列的情況，如果是m行1列就改成A(p)=[];
station(q)=[]; %對(duì)應(yīng)站點(diǎn)刪掉

%% 對(duì)缺測(cè)值進(jìn)行線性插值塌计，利用fillmissing函數(shù)
for i = 1:length(data_d2)
data_d2{1,i}(:,5) = fillmissing(data_d2{1,i}(:,5),'linear');
end
%如果是降雨等數(shù)據(jù)挺身，注意去掉負(fù)值
for i = 1:length(data_d2)
data_d2{1,i}(data_d2{1,i} < 0) = 0;
end

%% 單位轉(zhuǎn)換
for k = 1:length(data_d2)
temp = data_d2{1,k};
temp(:,5) = temp(:,5)/100; %本數(shù)據(jù)降水需除以10
data_d2{1,k} = temp;
end