衛(wèi)星導航定位算法及程序設(shè)計-實驗報告

1、-. z. - - . 可修編. 2013 級測繪工程專業(yè)衛(wèi)星導航定位算法與程序設(shè)計實驗報告實驗名稱：衛(wèi)星導航根本程序設(shè)計班級：*：*：實驗時間： 2016年6月28日2016年6月30 中國礦業(yè)大學目錄 TOC o 1-3 h z u HYPERLINK l _Toc455149467實驗一時空基準轉(zhuǎn)換 PAGEREF _Toc455149467 h 2HYPERLINK l _Toc455149468一、實驗目的 PAGEREF _Toc455149468 h 2HYPERLINK l _Toc455149469二、實驗內(nèi)容 PAGEREF _Toc455149469 h 2HYPERLI

2、NK l _Toc455149470三、實驗過程 PAGEREF _Toc455149470 h 2HYPERLINK l _Toc455149471四、實驗感想 PAGEREF _Toc455149471 h 2HYPERLINK l _Toc455149472實驗二 RINE*文件讀寫 PAGEREF _Toc455149472 h 2HYPERLINK l _Toc455149473一、實驗目的 PAGEREF _Toc455149473 h 2HYPERLINK l _Toc455149474二、實驗內(nèi)容 PAGEREF _Toc455149474 h 2HYPERLINK l _To

3、c455149475三、實驗過程 PAGEREF _Toc455149475 h 2HYPERLINK l _Toc455149476實驗三衛(wèi)星軌道計算 PAGEREF _Toc455149476 h 2HYPERLINK l _Toc455149477一、實驗目的 PAGEREF _Toc455149477 h2HYPERLINK l _Toc455149478二、實驗內(nèi)容 PAGEREF _Toc455149478 h 2HYPERLINK l _Toc455149479三、實驗過程 PAGEREF _Toc455149479 h 2HYPERLINK l _Toc455149480四、實

4、驗感想 PAGEREF _Toc455149480 h 2-. z.實驗一時空基準轉(zhuǎn)換一、實驗目的1、加深對時空系統(tǒng)及其之間轉(zhuǎn)換關(guān)系的理解2、掌握常用時空基準之間的轉(zhuǎn)換模型與軟件實現(xiàn)3、每人獨立完成實驗規(guī)定的內(nèi)容二、實驗內(nèi)容本實驗內(nèi)容包括：內(nèi)容一：編程實現(xiàn)GPS起點1980年1月6日0時對應(yīng)的儒略日內(nèi)容二：編程實現(xiàn)2011年11月27日對應(yīng)的GPS周數(shù)與一周內(nèi)的秒數(shù)內(nèi)容三：在WGS84橢球的條件下，編程實現(xiàn)當中央子午線為117度時，計算高斯坐標* = 3548910.811290287, y = 179854.6172135982 對應(yīng)的經(jīng)緯度坐標？內(nèi)容四：WGS84橢球下，外表*=-2408

5、000; y=4698000;z= 3566000處的地平坐標系坐標為: e=704.8615;n=114.8683;u=751.9771的點對應(yīng)的直角坐標為多少？三、實驗過程1.針對第一、二局部內(nèi)容：1.1解決思路：先建立 TimeStruct.h的頭文件，將格里高利歷、GPS時間構(gòu)造、儒略日時間構(gòu)造共構(gòu)造體的方式放在里面；在建立“TimeTr的頭文件，建立類“CTimeTr，創(chuàng)立變量“GPSTime、“Time、JulDay，并且申明函數(shù)“TIME2JUL、“TIME2GTIME等，用這些函數(shù)分別實現(xiàn)所需要的轉(zhuǎn)換。1.2具體的實現(xiàn)函數(shù)：“TIME2JUL函數(shù)：double CTimeTr:

6、TIME2JUL()/TIME Time,JULIANDAY &JulDaydouble m,y;double D;/h =Time.byHour+Time.byMinute/60.0+Time.dSecond/3600.00;if(Time.byMonth=2)y=Time.wYear-1;m=Time.byMonth+12;elsey=Time.wYear;m=Time.byMonth;D=floor(365.25*(y+4716)+floor(30.6001*(m+1)+Time.byDay+Time.byHour/24.0-1537.5;JulDay.lDay = int(D);Jul

7、Day.lSecond = D-int(JulDay.lDay);return 0;“TIME2GTIME：void CTimeTr:TIME2GTIME()double JD;long m,y;int WN;double Wsecend;/UT=Time.byHour+Time.byMinute/60.0+Time.dSecond/3600.00;if(Time.byMonth=0.0001 & abs(l1-l0)=0.0001 B0=B1; aa0=aa1; l0=l1; F_*B=(c*Beta2+(c*Beta4+(c*Beta6+c*Beta8*cos(B0)2)*cos(B0)2

8、)*cos(B0)2)*sin(B0)*cos(B0); F_*Bl=a2*l02+a4*l04+a6*l06; F_yBl=a3*l03+a5*l05; B1=(*-F_*B-F_*Bl)/(c*Beta0); aa1=(a*cos(B1)/sqrt(1-e2*(sin(B1)2); l1=(y-F_yBl)/aa1;endL=rad2deg(l1)+L0;B=rad2deg(B1);2.3實驗結(jié)果四、實驗感想本次試驗是花時間較多的一次實驗，關(guān)于時間轉(zhuǎn)換的局部全部都是自己動手將matlab代碼寫成“C+的類，進展實現(xiàn)的。其中遇到的較大的困難是儒略日向UTC轉(zhuǎn)換的局部，這局部的函數(shù)步驟較多，關(guān)

9、鍵是在一開場的時間構(gòu)造里面，各時間各局部的數(shù)據(jù)類型大多定義的是“int型的，但是在進展計算的時候有較多的小數(shù)，需要用到浮點型的函數(shù)，這局部用了較多的時間。在做這個實驗的時候，第一天花了時間主要是轉(zhuǎn)換代碼，使程序沒有錯誤，能夠正常的運行出來，出現(xiàn)黑框框，但是還只有個別功能能夠用，能夠運行出正確的結(jié)果；第二天時間主要是花在修改函數(shù)上頭，能夠使所寫的功能都能運行出正確的結(jié)果。通過做時間轉(zhuǎn)換的實驗，使自己產(chǎn)生了第一次親自編寫“C+代碼的經(jīng)歷，而且所有錯誤的解決全部都是自己解決，收獲不少。實驗二 RINE*文件讀寫一、實驗目的1、深入了解RINE*文件格式2、進一步提高MATLAB程序設(shè)計能力3、掌握N

10、文件、O文件、SP3文件的根本讀寫技巧二、實驗內(nèi)容本實驗內(nèi)容包括：1、任選IGS站，下載N文件、O文件與SP3文件；2、編程實現(xiàn)N文件讀入，并采用中文標注出主要參數(shù)的名稱及作用；編程實現(xiàn)O文件讀入，并采用中文標注出主要參數(shù)的名稱及作用；5、編程實現(xiàn)SP3文件讀入，并采用中文標注出主要參數(shù)的名稱及作用；三、實驗過程1、針對第一局部內(nèi)容：編程實現(xiàn)N文件讀入，并采用中文標注出主要參數(shù)的名稱及作用1.1、解決思路：按照“GPSeasy開源代碼提供的函數(shù)，按照實驗要求讀取了N文件的內(nèi)容，先用“rine*e函數(shù)，將N文件讀取成“eph.dat文件，然后再用“get_eph函數(shù)將“eph.dat文件讀取成“

11、Eph矩陣，此矩陣中包含了N文件中的數(shù)據(jù)，在最后用“fprintf函數(shù)將所需要的數(shù)據(jù)輸出成.T*T文件即可。1.2、主要函數(shù)代碼:“get_eph函數(shù)：function eph = get_eph(ephemeridesfile)fide = fopen(ephemeridesfile);eph, count = fread(fide, Inf, double);noeph = count/22;eph = reshape(eph, 22, noeph);“rine*e函數(shù)局部：function rine*e(ephemerisfile, outputfile)fide = fopen(eph

12、emerisfile);head_lines = 0;while 1 head_lines = head_lines+1; line = fgetl(fide); answer = findstr(line,END OF HEADER);if isempty(answer), break; end;end;head_lines主函數(shù)中輸出結(jié)果得函數(shù)局部：af0=data(19);%衛(wèi)星中差 M0=data(3); roota=data(4); deltan=data(5); ecc=data(6); omega=data(7); cuc=data(8); cus=data(9); crc=da

13、ta(10); crs=data(11); i0=data(12); idot=data(13); toe=data(18); af1=data(20);%對所要輸出的參數(shù)賦值fprintf(fid,n衛(wèi)星編號：%dn衛(wèi)星鐘差：%dn平近點角距：%dn軌道長半軸的平方根：%dn平均運動修正量：%dn軌道偏心率：%dn近地點角距：%dn緯度幅角的余弦調(diào)和項改正的振幅,prn,af0,M0,roota,deltan,ecc,omega,cuc);fprintf(fid,緯度幅角的正弦調(diào)和項改正的振幅：%dn軌道半徑的余弦調(diào)和項改正的振幅：%dn軌道半徑的正弦調(diào)和項改正的振幅：%dn軌道傾角：%dn

14、軌道傾角變化率：%dn星歷參考時刻：%dn,cus,crc,crs,i0,idot,toe)fclose(fid);1.3、輸出結(jié)果2、針對第二局部內(nèi)容：編程實現(xiàn)O文件讀入，并采用中文標注出主要參數(shù)的名稱及作用；2.1、實現(xiàn)思路：通過matlab的函數(shù)“fopen讀取O文件，得到O文件的指針，通過“anheader函數(shù)將文件中的接收機大致位置 appro*_*YZ1，天線的偏移值 ant_delta1，觀測值類型“Obs_types1等讀入成為matlab的矩陣，然后通過循環(huán)，利用“grabdata函數(shù)將所需要的歷元的觀測文件依次輸出來，最后通過“fprintf函數(shù)，將所需要的數(shù)據(jù)依次打印出來

15、。2.2、主要函數(shù)：“anheader函數(shù)：function Obs_types, ant_delta,ifound_types,appro*_*YZ = anheader(file)fid = fopen(file,rt);eof = 0;ifound_types = 0;Obs_types = ;ant_delta = ;appro*_*YZ = ;while 1 % Gobbling the header line = fgetl(fid); answer = findstr(line,END OF HEADER);if isempty(answer), break; end;if (l

16、ine = -1), eof = 1; break; end; answer = findstr(line,ANTENNA: DELTA H/E/N);if isempty(answer)for k = 1:3 delta, line = strtok(line); del = str2num(delta); ant_delta = ant_delta del;end;end answer = findstr(line,APPRO* POSITION *YZ);if isempty(answer)for k = 1:3 app_*YZ, line = strtok(line); del = s

17、tr2num(app_*YZ); appro*_*YZ = appro*_*YZ del;end;end answer = findstr(line,# / TYPES OF OBSERV);if isempty(answer) NObs, line = strtok(line); NoObs = str2num(NObs);for k = 1:NoObs ot, line = strtok(line); Obs_types = Obs_types ot;end; ifound_types = 1;end;end;%fclose(fid);“grabdata函數(shù)：function Obs =

18、grabdata(fid, NoSv, NoObs)%GRABDATA Positioned in a RINE* file at a selected epoch% reads observations of NoSv satellitesglobal linObs = zeros(NoSv, NoObs);if NoObs = 5 % This will typical be Turbo SII datafor u = 1:NoSv lin = fgetl(fid);for k = 1:NoObs Obs(u,k) = str2num(lin(2+16*(k-1):16*k-2); end

19、endelse% This will typical be Z12 data Obs = Obs(:,1 2 3 4 5); % We cancel the last two columns 6 and 7 NoObs = 5;for u = 1:NoSv lin = fgetl(fid); lin_doppler = fgetl(fid);for k = 1:NoObs %-1if isempty(str2num(lin(1+16*(k-1):16*k-2) = 1, Obs(u,k) = nan; else% Obs(u,k) = str2num(lin(1+16*(k-1):16*k-2

20、); end% Obs(u,NoObs) = str2num(lin(65:78);endendend2.3實驗結(jié)果四、實驗感想這局部實驗是我在之前做的，之前自己有看過“gps_easy有關(guān)的代碼，看過相關(guān)的“N文件“O文件讀寫函數(shù)，并且學會了如何調(diào)用這些函數(shù)，對里面的輸出量有了一點的了解，所以我自己的主要工作就是運用了“fprintf函數(shù)，將讀取到matlab中的矩陣寫入T*T文檔中，這局部工作量不是很大，但較有意義。實驗三衛(wèi)星軌道計算一、實驗目的1、進一步熟悉N文件的讀入2、掌握開普勒參數(shù)計算衛(wèi)星軌道的過程3、編程實現(xiàn)采用播送星歷計算衛(wèi)星軌道4、掌握MATLAB函數(shù)調(diào)用步驟二、實驗內(nèi)容本實

21、驗內(nèi)容包括：1、調(diào)試時間轉(zhuǎn)換函數(shù)，熟悉內(nèi)容，備主函數(shù)調(diào)用2、調(diào)試播送星歷導航文件的讀入程序，備主函數(shù)調(diào)用3、根據(jù)衛(wèi)星位置計算公式編寫主函數(shù)，同時調(diào)用時間轉(zhuǎn)換、星歷讀取等的子函數(shù)來共同完成衛(wèi)星位置的計算，最后輸出計算結(jié)果4、理清程序各模塊的功能構(gòu)造三、實驗過程1、實驗思路：在教師提供的“SPP文件中，直接利用衛(wèi)星位置計算函數(shù)，進展衛(wèi)星位置的計算，將利用“Gps.cpp文件中的 GetGpsPosition函數(shù)，利用其中的迭代求解衛(wèi)星位置局部，用“cout直接將衛(wèi)星迭代后的位置直接輸出，因為星歷文件中有較多的星歷，所以利用循環(huán)語句，將求解出來的衛(wèi)星位置依次輸出出來。2、主要函數(shù)bool CGps:

22、GetGpsPosition()GPSTIME ts;GPSTIME tr;GPSTIME ts0;GPSTIME oGTime;GPSTIME nGTOC;vectorSendSignPosition;GpsPos GpsPTemp;GpsSendPosition SdSignPoTemp;int nTheFitPoint=0;Position pTemp;GpsPosition.clear();if (oData.size()=0)strErr = (PRN=%d沒有對應(yīng)星歷);returnfalse;if (nData.size()=0)returnfalse;coutfi*ed;for (int i=0;ioData.size();i+)double *V*=newdouble4;/for (int j=0;j1e-007);/earthrot(tr.lSecond-ts.lSecond,pTemp.*,pTemp.YY,pTemp.ZZ);TimeToGpsTime(nDatanTheFitPoint.TOC,nGTOC);/SdSignPoTemp.dt=nDatanTheFitPoint.d