基于MATLAB的标准化降水指数SPI计算程序

%计算各站点SPI1，SPI3SPI6SPI9SPI12及季SPI、年SPI值
%改变line7—读取excel数据（只包括站点和日期包围的降雨数据）
%若计算其他尺度如SPI24，改变line9 （数据尺度） ； 改变line66-68
%改变line74—78（输入到excel中）
%code by cc 2018-9-2
close all;clc;
rainData = xlsread（‘E:\weihe‘1‘B5:P640‘）; %读取降雨数据
[MN] = size（rainData）;  %m总月数n站点数
years=M/12;  scale=[136912];
[spi1spi3spi6spi9spi12]=deal（zeros（MN））; %136912各站点的spi值    
[spiSprspiSumspiAutspiYear]=deal（zeros（yearsN））; spiWin=deal（zeros（years-1N））; %四季各站点的spi值
time=0;
%%%%%%%%%%%%%%fList=zeros（M1）; mList=zeros（M1）; nList=zeros（M1）;%%%%%%%%%
for station=1:N
    rainDataSta=rainData（:station）;       %某一站降雨数据
    spiArrayT=zeros（Mlength（scale））;      %某站5个尺度spi矩阵
    sheet=1;
    for timeScale=scale                    %分不同时间尺度
        rainAdd=zeros（M-timeScale+11）;    %时间尺度下的累计雨量序列
        num=length（rainAdd）;
        spiColumn=zeros（M1）; 
        sheet=sheet+1;
        for i=timeScale:M
            for j=0:timeScale-1
                rainAdd（i-timeScale+1）=rainAdd（i-timeScale+1）+rainDataSta（i-j）;
            end
        end 
        %#################################################
        %计算12个序列下SPI值所需的参数和降雨为0的概率
        qList=zeros（121）;para1=zeros（121）;para2=zeros（121）;
        for month=1:12     %列年里12个同期序列
            rainSameT=rainAdd（month:12:num）;   %某时间尺度下同期降雨序列
            number=length（rainSameT）;
            %计算降雨量为0的概率
            m=0 ; n=0;
            for i=1:number
                if rainSameT（i）==0
                    m=m+1;
                else
                    n=n+1;
                end
            end
            q=m/number;
            rainNot0=rainSameT（rainSameT~=0）;