MATLAB代做太阳黑子模拟分析.docx

《MATLAB代做太阳黑子模拟分析.docx》由会员分享，可在线阅读，更多相关《MATLAB代做太阳黑子模拟分析.docx（9页珍藏版）》请在第一文库网上搜索。

1、太阳黑子是人们最早发现也是人们最熟悉的一种太阳表面活动。因为太阳内部磁场发生变化，太阳黑子的数量并不是固定的，它会随着时间的变化而上下波动，每隔一定时间会达到一个最高点，这段时间就被称之为一个太阳黑子周期。太阳黑子的活动呈现周期性变化是由施瓦贝首次发现的。沃尔夫（RWo1fer）继而推算出11年的周期规律。实际上，太阳黑子的活动不仅呈11年的周期变化，还有海耳在研究太阳黑子磁场分布时发现的22年周期；格莱斯堡等人发现的80年周期以及蒙德极小期等。由于太阳黑子的活动规律极其复杂，时至今日科学家们仍在努力研究其内在的规律和特性。事实上，对太阳黑子活动规律的研究不仅具有理论意义,而且具有直接的应用需

2、求。太阳黑子的活动呈现周期性变化的，沃尔夫（RWo1fer）根据在过去的288年（1700年7987年）间每年太阳黑子出现的数量和大小的观测数据推算出11年的周期规律。我们利用MatIab强大的数据处理与仿真功能，对Wo1fer数进行功率谱密度分析从而可以得到对太阳黑子活动周期的结论。Basicproject：。Identifyareasonab1efrequency,题目解析：这个题目的要求就是计算出太阳黑子的大致的活动频率。解决方法：导入数据，得到如下的波形：SSN这个步骤计算得到的太阳黑子的运动周期为:太阳黑子的活动周期为(单位年):JnainPeriod=11.1339这个步骤的原理，

3、是通过简单的FFT估计频谱图，得到最后的频率值，从而计算周期。通过计算FFT得到频谱图，然后计算频谱的最大值，其倒数为周期。ComputeX1=cos(2ft)and&=sin(2).题目解析：这个题目的含义就是根据第一步计算得到的频率来计算X1和X2。解决方法：直接代入f,X1=CoS(2；FmX2=Sin(2;Tm在MAT1AB中，其代码如下所示:40-t=1/(12):1/(12):Iength(SSN)/(12);41 -yc=cos(2*pi*tmainPeriod);42 -ys=sin(2*pi*tainPeriod);43 -figure;44 -subp1ot(R145 -P

4、IQtiyeJnh1f，146 -p1ot;jho1doff47 -IegendCX1,U);,48 -49 -p1ot(yc,r一,);ho1don50 -p1ot(ys,k-);ho1don51 -p1ot(SSN/max(SSN),b.,);ho1doff52 -1egend(,X1,X2,SSN,)：其仿真结果如下所示:。Fitthe1inearizedcyc1icmode1with1eastsquaresusingMat1ab.题目解析：使用最小二乘法得到线性方程模型的估计结果解决方法：这里，使用的线性模型为：N,=A)+B1X1+A+error这个部分的MAT1AB代码为:61一6

5、2-63-64-65一66一67-68-69-70一71一72一737475-76-77fork=1:Iength(SSN);%开始求Kh1=1,yc(k),ys(k)：X=h*p*h1+99：x1=inv(x):k1=p*h1*x.1;d1=SSN(k)-h1,*c;c1=c+k1*d1;%开始求%求出KF%求被辨:hS1圜K(k)向值识参数C当前值与上一次的值的差值的相对变化相对变化的列向量加入误差矩阵的最后一列e1=c1-c;e2=e1./c;e(:,k)=e2;Hh()gic%把当前c=111;%直接给出被辨识参数的初始值，即一个充分小的实向量p=10A6*eye(X3);%直接给出初

6、始状态PO,即一个充分大的实数单位矩阵c=c1;:然蜃号的叁数作为下一次递推的旧参数c(:,k)=c1;靠我参数C列向量加入辨识参数矩阵的最后一列p1=pO-k1*k*(h*p*h1+1;%求出p(k)的值p=p1;%给下次用end代码中CO(1),C0(2),CO(3)分别为估计得到的系数/),4,尸2将上面的图放大之后，得到的局部估计效果如下图所示:OInterpretthefittedmode1parametersintermsofamp1itude,frequency,andphaseshiftofthecyc1icmode1.题目解析：根据上面的最小二乘法得到的结果计算模型的幅度，频

7、率和相位。解决方向：t=BCOS(2万ft+)=Pcos(2万力)cos()一sin(2ft)sin(/)=COS(夕)COS(2万ft)-sin(%r)sin(2ft)=iX+2X1=arctan(-y2)=arctan(夕Si11()夕=JQ+4这个部分的代码如下所示:94-95-96-97-98-99-IOO-101102103104-105-t=1/(12):1/(12):1ength(SSN)Z(12)；fork=1:Iength(SSN);beta(k)beta1(k)beta2(k)BETA(k)fai(k)=c(1,k)；=c(2,k)；=ck);=sqrtJ=atanV%ut

8、=Bcos(2*pf*-jJ;%t=%yc=Ut(k)1/(12):1/cos(2*p*tmaTnPerod);=BETA(k)*cos(2*pi*t(k)mainPeriod+fai(k)+betaO(k);end通过上面的步骤，得到4=/?cos（2#+夕）中的几个参数，其仿真结果如下所示:QUseresidua1p1otstova1idate/eva1uatethemode1.题目解析：这个步骤就是画出残差图：解决方向：Considerhowthemode1mightbeimproved.题目解析：模型的改进解决方向：这里，主要引入多个变量作为参数估计，这里，将最小二乘法的估计函数为：M

9、=A)+t+2cos(24ft)+3sin(2乃ft)+ycos(4Z)+5CoS(4万ft+sin(6%ft)+1COS(6万力)由于从历年的数据可知，虽然总体上数据呈现出周期波形，但是对于每个周期中的波形,其中的局部图中会出现部分高频周期量，所以在估计的时候，这里加上更高频率的分量。使用和上面所讲的最小二乘法进行上述参数的估计。4505152554555657fork=1:1ength(SSN);%开始求Ke(:,k)c1,k,yc(k),ys(k),yc2(k),ys2(k),yc3(k),ys3(k):h*ph1+99:inv(x):p*h1*x1;=SSN(k)-h1,*cO;=c+

10、k1dkNJjJc1-c.=e1,c.=。2:c(:,k)c1%开始求K(k)X求出K的值ighSD堡寺通前值与上一次的值的差值-THr求学敷的相对变化OiI0汜当前相对变化的列向量加入误差矩阵的最后一列Q笔静播的参数作为下一次递推的旧参数QQ12248W躺忌数C列向置加入辨识参数矩阵的最后一列p1=pO-k1*k*(h*pO*h1+1;X求出p(k)的值p=p1;福下次用end%估计得到的太阳黑子活动曲线fork=1:1enth(SSK);X开始求KY_predict2(k)c(1,k)+c(2,k)k+c(3,k)*yc(k)+c(4,k)ys(k)+c(5,k)*yc2(kend这个步骤，得到的仿真结果如下图所示: