太阳黑子活动周期的分析

引言太阳黑子是人们最早发现也是人们最熟悉的一种太阳表面活动。因为太阳内部磁场发生变化，太阳黑子的数量并不是固定的，它会随着时间的变化而上下波动，每隔一定时间会达到一个最高点，这段时间就被称之为一个太阳黑子周期。太阳黑子的活动呈现周期性变化是由施瓦贝首次发现的。沃尔夫(R.Wolfer)继而推算出11年的周期规律。实际上，太阳黑子的活动不仅呈11年的周期变化，还有海耳在研究太阳黑子磁场分布时发现的22年周期；格莱斯堡等人发现的80年周期以及蒙德极小期等。由于太阳黑子的活动规律极其复杂，时至今日科学家们仍在努力研究其内在的规律和特性。事实上，对太阳黑子活动规律的研究不仅具有理论意义，而且具有直接的应用需求。太阳黑子的活动呈现周期性变化的，沃尔夫(R.Wolfer)根据在过去的288年(1700年~1987年)间每年太阳黑子出现的数量和大小的观测数据推算出11年的周期规律。我们利用Matlab强大的数据处理与仿真功能，对Wolfer数进行功率谱密度分析从而可以得到对太阳黑子活动周期的结论。实验原理在该实验中，对Wolfer数序列做FFT变换后得到Y（长度为n），只取其前个数据的功率谱密度的估计值。原因是时域为离散的实序列的傅立叶变换对应于具有周期性且偶对称的频域特性，因此Y的前个数据已经包含了Wolfer数的全部信息。根据DFT的频域单位与DTFT的频域单位的表达关系式以及与f对应关系，可以看出与f呈线性关系；同样地，因为Y的前个数据已经包含了Wolfer数的全部信息，只取前个数据分析功率-频率图时，对应的横坐标时应取。（注：为采样频率）三、方案组成与流程图111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111四、实现步骤及内容(1)读取数据并绘制Wolfer图，观察太阳黑子的活动现象。程序如下：loadsunspot_sun.dat%读取数据year=sunspot(:,1);%读取年份信息wolfer=sunspot(:,2);%读取黑子活动数据figure%新建图像plot(year,wolfer)%画出时域图xlabel('Years');ylabel('SunspotData');title('SunspotData')%标注横纵坐标和标题pause为清楚起见，画出部分时间段(如实验中采取开始50年)的Wolfer图进行观察。程序如下：figureplot(year(1:50),wolfer(1:50),'b.-');%用前50年数据做图xlabel('Years');ylabel('SunspotData');title('Atthefirst50years')%标注横纵坐标和标题(2)对已经得到的Wolfer数应用FFT技术分析它的变化规律。程序如下：Y=fft(wolfer);%对全部数据做FFTY(1)=[];%由于算法问题，Y(1)为所有fft数值之和，因此需舍弃第一个点pause观察Wolfer数在FFT分析后得到的在复平面中的分布图。程序如下：figureplot(Y,'ro')%在复平面做图，空心点title('FourierCoefficientsintheComplexPlane');xlabel('RealAxis');%标注横坐标为’实部’ylabel('ImaginaryAxis');%标注纵坐标为’虚部’pause(3)功率与频率的关系曲线则被定义为周期图。用周期图法对Wolfer数进行功率谱密度分析。程序如下：figuren=length(Y);%取FFT结果长度为npower=abs(Y(1:n/2)).^2;%取前n/2个数据求其功率nyquist=1/2;%取最大频率为0.5freq=(1:n/2)/(n/2)*nyquist;%将FFT转换为DTFTplot(freq,power)%画周期图xlabel('cycles/year');title('Periodogram')pause为清楚起见，对功率和频率的前50个分量作它的周期图。程序如下：figureplot(freq(1:50),power(1:50))%为观察方便，取前50年数据做图xlabel('cycles/year')pause(4)确定出太阳黑子的活动周期。为清楚起见，画出功率与周期（频率的倒数）的关系曲线图。程序如下：figureperiod=1./freq;%将横坐标单位转换为年/周期plot(period,power);%画周期图axis([05002e+7]);ylabel('Power');xlabel('Period(Years/Cycle)');%单位为年/周期pause(5)在功率与周期关系曲线图中标出功率的最高点，该位置对应的周期则为太阳黑子活动周期。holdon;index=find(power==max(power));%找到频率最大点，该点横坐标即为太阳黑子周期mainPeriodStr=num2str(period(index));plot(period(index),power(index),'r.','MarkerSize',25);%用实心点指出该点text(period(index)+2,power(index),['Perio