1、AssignmentofAssignmentofInertialTechnology慣性技術(shù)作業(yè)（2016秋）MyChineseNameMyStudentNo.16S104Thereportistocontain:1. Descriptionofthetasks-contentsofthenexttwopagesandthepreviousassignments.2. Thecodeofyourprograms,andtheirexplanation.3. Theresultsofyourcomputationorsimulation(aslistedbytherequirement).4.

2、Youranalysisoftheresult,andyourreflectionontheprogrammingorsimulation5. Originalitystatementsorreference/assistanceacknowledgements.Englishisexpectedinwriting,thoughChineseisalsoaccepted.Assignment1:2-DOFresponsesimulationAssignment1:2-DOFresponsesimulationI.DescriptionofthetasksI.Descriptionoftheta

3、sksA2-DOFgyrohasarotorwithangularmoment10000gcms.Itsequatorialangularinertiasareboth4gcms2.Pleaseinvestigatetheresponseofthegyrotothefollowingtypesoftorquesaslistedinthetable,andpresentwhateveryoucandiscoverorconfirmfromtheresults.Table1-1TypesofinputtorquesTorquedescriptionDirectionSuggestedsimulat

4、iontimeImpulse,magnitude40000gcm,duration1e-5sInnerringaxis0.0030.02sConstant,magnitude1gcm.Innerringaxis0.0050.02sSinusoidal,amplitude1gcm,frequency2HzOuterringaxis0.52sSinusoidal,amplitude1gcm,frequency5HzOuterringaxis0.20.8sSinusoidal,amplitude1gcm,frequency10HzOuterringaxis0.10.4sSinusoidal,ampl

5、itude1gcm,frequency20HzOuterringaxis0.050.4sSinusoidal,amplitude1gcm,frequency50HzOuterringaxis0.020.2sNotethatthedefaultconfigurationsofthesimulinkparameters,suchasthoseofstepsize,mightnotbeadequatetobringaboutthenuancesinthesimulationresult.2.SimulationandAnalysis2.SimulationandAnalysisAccordingto

6、thetransferfunctionofthe2-DOFgyro,theoutputscanbeexpressedas:Giventhisfunction,wecanestablishtheblockdiagramofthesysteminSimulink.2.1Input1:Impulse,magnitude40000g.cm,duration1es,oninnerringaxis2.1Input1:Impulse,magnitude40000g.cm,duration1es,oninnerringaxisInthiscase:My=40000-(t),Mx1=0Thefrequencyo

7、fthenutationisobtainedas:Radiusofnutationisobtainedas:My*410rad0.1375arcminHTheblockdiagramofthesysteminSimulinkisshowedinFig2.1.JyJxJys2H2MXIXIs-sJxJys2H2MysJxJxJys2H2MyssJxJys2H2Mx1sH100002500rad/s397HzJe內(nèi)環(huán)肽沖響應(yīng)Fig2.1TheblockdiagramofinputlinSimulinkIntheblockdiagram,wecanaddaXYGraphtoshowtherelati

8、onshipnbetWeenowever,thegraphoftheirrelationshipcannotbeedited.SoIlogthedataofScopeandScope1intoworkspaceofMatlabandnamethemseparatelyasAlphaandBeta.Firstly,Isetthesimulationstoptimeas0.003s,thenrunthesystem.Inworkspace,IdrawthegraphofthetworesponseintooneaxesandgetthefollowingFig2.2.ThenIalsodrawth

9、etrajectoryof2-DOFgyrosresponsetothisimpulseinputinFig2.3.Fig2.22-DOFgyrostworesponsestoimpulseinput(0.003s)FromFig2.4wecanseethatthetrajectoryisnotsmooth,why?Foronethingwecanconfirmisthatitisrelatedtothesimulationtime.BecausetheonlydifferencebetweenFig2.4andFig2.3isthesimulationtime.Atthebeginning,

10、wehavefiguredoutthatJe100002500rad/s397Hz4Whichmeansthattheperiodis1T=-=0.0025s:0.003sf397So0.003sisaboutoneperiodand0.02sisaboutnineperiods.OwingtothesimulationstepofSimulinkisfixed(ormaybeIdontknowhowtochange),theresultwillhasvwIbensimulationtimeislong,whichmayexplainthesawtootharoundtheedge.2.2In

11、put2:Constant,magnitude1g.cm,oninnerringaxis2.2Input2:Constant,magnitude1g.cm,oninnerringaxisInthiscase:My=1回(t),M*=0Thefrequencyofthenutationisobtainedas:H10000=2500rad/s397HzJeThedriftrateisobtainedas:My4=1104rad/s0.345deg/in=20.7deg/Theradiusofnutationisobtainedas:JeMyHT=410-8rad1.3610“arcsecTheb

12、lockdiagramofthesysteminSimulinkisshowedinFig2.5.ThetworesponsestotheconstantinputFig2.3Trajectoryof2-DOFgyrosresponsetoimpulseinput(0.003s)FromFig2.3wecanobviouslygetthatthetrajectoryof2-DOFgyrosresponsetoimpulseinputisacircle,whosecenteris(-410-,0)rad,andradiusis410-5rad.Inordertocomparethedatawit

13、hdifferentsimulationtime,Ichangethesimulationstoptimeto0.02s,andgetthetrajectoryinFig2.4.inTrujectar/of2-DOFgroiresponsetuimpulseinputAtphrad)Atphrad)Fig2.4Trajectoryof2-DOFgyrosresponsetoimpulseinput(0.02s)within0.01sisshowedinFig2.6.Thetrajectoryof2-DOFgyrosresponsetothisconstantinputisshowedinFig

14、2.7.內(nèi)環(huán)階躍響應(yīng)Fig2.5Theblockdiagramofinput2inSimulink(PBSflJJBACJ5C0(PBSflJJBACJ5C0Fig2.62-DOFgyrostworesponsestoconstantinput(0.01s)Fig2.7Trajectoryof2-DOFgyrosresponsetoconstantinput(0.01s)FromFig2.7wecanobviouslygetthatthetrajectoryof2-DOFgyrosresponsetoconstantinputisacycloid.2.32.3 Input3:Sinusoida

15、l,amplitudeIg.cm,frequency2Hz,onouterringaxisInput3:Sinusoidal,amplitudeIg.cm,frequency2Hz,onouterringaxisInthiscase:Mx1=16in(4二t),My=0Thefrequencyofthenutationisobtainedas:TheblockdiagramofthesysteminSimulinkisshowedinFig2.8.Thesimulationresultsintimedomainareshowedinthefollowingfigures.Fig2.9isthe

16、outputof2-DOFgyrostworesponseswithin2s.Fig2.10istheoutputofouterringwithin2s.Fig2.11isthetrajectoryof2-DOFgyrosresponsetosinusoidalinputwithin0.5s.AswecanseefromtheFig2.10,thereareobvioussawtoothwaveintheoutputoftheinnerring.Itunexpectedphenomenoninmyoriginaltheoreticalanalysis.0Je10000:2500rad/s:39

17、7Hz外環(huán)正弦響應(yīng)正獨(dú)刀犯內(nèi)印信函P PFig2.8Theblockdiagramofinput3inSimulink(p何)口3 3一加A0如cQdmfucQdmfuH HE E.W.WFig2.10Theoutputofouterringfromsinusoidalinput(2s)Fig2.11Trajectoryof2-DOFgyrosresponsesinusoidalinput(0.5s)Ibelievethesawtoothwaveiscausedbythenutation.Forthefrequencyofthenutationisobtainedas:whichisfarhi

18、gherthanthefrequencyoftheappliedsinusoidaltorque,namely:wa AlphB(rad)AlphB(rad)Fig2.19TheoutputofouterringfromFig2.20Trajectoryof2-DOFgyrossinusoidalinput(0.02s)responsetosinusoidalinput(0.02s)FromFig2.8wecangettheoutputofouterringintimedomainasshowedinFig2.19.Fig2.17TheoutputofouterringfromFig2.18T

19、rajectoryof2-DOFgyrossinusoidalinput(0.05s)responsetosinusoidalinput(0.05s)AJOhfloutputx xlOlOL,L,TrAjMTrAjM：(0fvd24X)FgymrMporiM(0fvd24X)FgymrMporiMtotosinusoidalinpiMsinusoidalinpiM3.Conclusion3.Conclusiontrajecteayohat-tDOFgyrosresporiseaprecisecirclewithanimpulseinputtorque,andcycloidwithaconsta

20、ntinputtorque.However,thecirclebecomesellipseandwithsawtoothwhentheinputtorquebecomesasinusoidalsignal.Besides,whenfrequencyofsignalincreases,thetrajectorybecomesmoreirregular,thatmaybebecausetheouterringsvibrationincludesAftercomparingthesimulationresults,itmorenoise.Assignment2:Single-axisINSsimul

21、ationAssignment2:Single-axisINSsimulationI.DescriptionofthetasksI.DescriptionofthetasksInafictitioustestofamagneticlevitationtrainalongatrackrunningnorth-south,itfirstacceleratesandthencruisesataconstantspeed.Onboardisasingle-axisplatformINS,workinginthewaydescribedbythecoursewareofUnit5:Basicproble

22、msofINS.ThemotioninformationandEarthparametersareshownintable2-1,andthepossibleerrorsourcesarelistedinTable2-2.YouareaskedtosimulatetheoperationoftheINSwithin9,900seconds,andinvestigate,firstonebyoneandthenaltogether,theimpactoftheseerrorsourcesontheperformanceoftheINS.Notethattheblockdiagraminthele

23、cturenotes(figure3.6ofboth2015and2016versions)ortheoldversionsofcoursewarehastobeslightlymodifiedbeforeyoucanobtainreasonableresults.Table2-1MotioninformationandEarthparametersMotioninformationvaluesunitsEarthparametersvaluesunitsInitialvelocity,northward5m/sAccelerationofgravity9.82m/sInitialpositi

24、on0mRadiusoftheEarth6371kmAcceleration,fromstart2,2m/sDurationofacceleration80sTable2-2PossibleerrorsourcesTypesvaluesunitsTypesvaluesunitsInitialpositionerror20mAccelerometerscalefactorerror0.00051Initialvelocityerror0.05m/sGyroscopescalefactorerror0.00051Initialplatformmisalignmenterror1Gyroscoped

25、riftingerror0.01o/hAccelerometerbiaserror0.000022m/s2.SimulationandAnalysis2.SimulationandAnalysisThereisonecorerelevantformula,togetthespecificformofitssolution,weshouldsubstitutetheunknownparameters.y.c=AN(1Ka)(y-ga)Firstly,theinputsignalisaccelerometeroftheplatform,andthevelocityoftheplatformisth

26、eintegrationoftheacceleration.p=-y/RTheaccelerationalongYpmaycontainstwoparts:fyp=ycos：-gsin二y-g：Whenaccelerometererrorsareconcerned,theoutputofaccelerometerwillbe:aN=(1Ka)fypAWhengyroerrorsconcerned:Only：isunknown:tt(1Kg)p；dt-；t0t全部誤差源Fig2.2TheSimulinkblockdiagramforAssignment2withallerrorsources.一

27、，pdtAcc.Av回目毋案wgyroi/w1十K.s splatformcomputercomputerFig2.1Thereferenceblockdiagraminthecourseware(rectified)kli&gratjrSkli&gratjrS群噓灣祐蜂IrivurabrZIrivurabrZ此此M Mo真五位置計(jì)策應(yīng)提AddAdd2二JoJo3PI I醯度十喀i i踏1小曰kiteqr加ni-算加速度算加速度聯(lián)敏lii&grdlDD3000013002(HD3000.歐S00090007O0UJUWHOBS00090007O0UJUWHOBFig2.

28、8Positionbiasoutputwhenonlyinitialplatformmisalignmenterrorexists0IOOO2tKKWOO500C60MMOUHUW加UUFig2.9PositionbiaswhenonlyaccelerometerbiaserrorexistsFig2.10PositionbiaswhenonlyaccelerometerscalefactorerrorexistsFig2.11PositionbiasoutputwhenonlygyroscalefactorerrorexistsFig2.12positionbiasoutputbiaswhe

29、nonlygyrodriftingerrorFig2.13PositionbiasoutputconsideringallerrorsourcesFig2.14PositionbiasoutputconsideringnoerrorsourcesAswecanseeintheabovesimulationresults,ifthereisnoerrorwecannavigatethetrainsmotioncorrectly,whichcomesfromnorthtothesouthasshowninFig2.3toFig2.5,beginningwithaconstantaccelerati

30、onwithin80secondsthencruisesataconstantspeed,approximately165m/s.However,thesituationwillchangealotwhendifferenterrorsputintothesimulation.TheinitialpositionerrorAy0=20meffectsleastasFig.2.6,forthiserrordoesntenterintotheclosedloopanditwontinfluencetheiterativeprocess.Thepositionbiasisconstantandcan

31、benegligible.Inthesecondcase,whentheaccelerometerscalefactorerrorexists,Ka=0.0005,asshowninFig2.10,theresultarestableandalmostaccurate,thepositionbiasisasinusoidaloutput.Soitiswiththe2accelerometerbiaserrorsituation,.：AN=0.00002m/s,inFig2.9,theinitialvelocityerror,2.y0=0.05m/sinFig2.7,andtheinitialp

32、latformmisalignmentangle,二0二1,inFig2.11.However,theinfluencedegreesofthedifferentfactorsarenotinthesamemagnitude.Theaccelerometerscalefactorinfluencestheleastwithmagnitudeof62,thentheinitialvelocitysmallermagnitudeof40,andtheaccelerometerbiasmagnitudeof25.Theinfluenceoftheinitialplatformmisalignment

33、angleismuchmoresignificantwithamagnitudeof62.Allthenavigationbiasinthesecondkindcaseissinusoidal,whichmeanstheyrelimitedandnegligibleastimepassesby.Inthethirdcase,suchasthegyroscalefactorerrorsituation,Kn=0.0005,inFig2.11,andthegyrogdriftingerror,；=0.01/h,resultsinFig2.12,effectsthemostsignificant,t

34、hetrajectoryofthenavigationdeviationaccumulatedastimegoes.Thepositionbiasisacombinationofsinusoidalsignalandrampsignal.Theyalsoshowthatthelongitudinalanddistanceerrorsresultedfromgyrodriftsarenotconvergentintime.Itmeanstheerrorsinthegyroscopedomostharmtoournavigation.Andduetothesignificantinfluenceo

35、fthegyrodriftingerrorsandthegyroscopescalefactorerror,resultsconsideringalltheerrorsources,andthenavigationpositionofthemotionwillbeawayfromtherealmotionafteranenoughlongtime,asshowninFig2.12.Thegyrodriftingerroristhemostsignificanteffectfactorofallerrors.Bythetimeof9900s,ithasreaches2250m,anditsnea

36、rlythequantityofthepositionbiasconsideringallerrorsources.Throughcontrastingalltheresults,wecanconcludethatthegyrodriftingerroristhemaincomponentofthewholepositionbias.3.Conclusion3.ConclusionThroughcontrastingalltheresults,wecanconcludethatthegyrodriftingerroristhemaincomponentofthewholepositionbia

37、s,andthegyrobiasorthedrifterrordomostharmtoournavigation.Soitisamustforustoweakenoreliminateitanyway.Inspiteofallthedisadvantagesdiscussedabove,theINSstillshowsusarelativelyaccurateresultsofsingle-axisnavigation.Assignment3Assignment3I.DescriptionofthetasksI.DescriptionofthetasksInanfictitiousmissio

38、n,aspaceshipistobeliftedfromalaunchingsitelocatedat19037NLand110o57EL,intoacircularorbit400kilometershighalongtheequator.ThespaceshipisequippedwithastrapdownINSwhosethreegyros,GX,GY,GZ,andthreeaccelerometers,AX,AY,AZ,areinstalledrespectivelyalongtheaxesXb,Yb,Zbofthebodyframe.CaseCase 1:1:Stationaryt

39、estStationarytestDuringapre-launchinggroundtest,thebodyframeofthespaceshipinitiallycoincideswiththelocalgeographicalframe,withitspitchingaxisXbpointingtotheeast,rollingaxisYbtothenorth,andheadingaxisZbupwards.Thenthebodyofthespaceshipismadetorotatein3steps:(1)80oaroundXb(2)90oaroundYb(3)1700aroundZb

40、Afterthat,thebodyofthespaceshipstopsrotating.Youarerequiredtocomputethefinaloutputsofthethreeaccelerometersinthespaceship,usingquaternionandignoringthedeviceerrors.Itisassumedthatthemagnitudeofgravityaccelerationatthegroundlevelisg0=9.79m/s2.CaseCase 2:2:ThelaunchingprocessThelaunchingprocessThespac

41、eshipisinstalledonthetopofaverticallyerectedrocket.Itsinitialheading,pitchingandrollingangleswithrespecttothelocalgeographicalframeare-90,90and0degreesrespectively.Thedefaultrotation|sequenceisheadingfpitchingfrolling.Thetoltherocketisinitially100mabovethesealevel.Thentherocketisfiredup.Theoutputsof

42、thegyrosandaccelerometersinthespaceshiparebothpulsenumbers.Eachgyropulseisanangularincrementof0.01arcsec,andeachaccelerometerpulseis1e-7g0,withg0=9.79m/s2.Thegyrooutputfre-quencyis100Hz,andtheaccelerometeris5Hz.Theoutputsofthegyrosandaccelerometerswithin1800sarestoredinaMATLABdatafilenamedmission.ma

43、t,mission.mat,containingmatricesGGMGGMof1800003fromgyrosandAAMAAMof90003fromaccelerometersrespectively.TheformatofthedatainthetwomatricesisasshowninthetabTheEarthcanbeseenasanidealsphere,withradiusR=6371.00kmandspinningratew=7.29210-rad/s,Theerrorsofthegyrosandaccelerometerscanbeignored,buttheeffect

44、ofheightonthemagnitudeofgravityhastobetakenintoaccount.Thegravityaccelerationatthesealevelofthenearequatorregioncanbechosenasg0=9.79m/s2,andthemagnitudeofgravityaccelerationataheightofhcanbecomputedasGXGYGZ-106584145-106487147-106686148-106486147-106687146-106585145-106587146-106585145-106784147-106

45、587147AXAYAZ-159846814634887-1279082-161193414636446-1253706-162539314637974-1228346-163885614639472-1202978-165231914640940-1177617-167562014635359-1221229-168684614636804-1196141-169807114638220-1171047-170930114639605-1145960-172054414640958-1120882les,with10rowsofeachmatrixselected.Eachrowrepres

46、entstheoutputsofthetypeofsensorsatasamplingtime.J?-Besides,theinfluenceofheightontheangularratesofthegeographicalframeandthechangingratesoflatitudeandlongitudeshouldalsobeconsidered.Velocity,positionandthegeographicalframecanbeupdatedevery0.2s,withinwhichtheattitudeofthevehiclecanbeupdatedmultipleti

47、mes,dependingonthechosenalgorithm(20for1-S,10for2-S,and5for4-S).Youarerequiredto:(1) computethefinalattitudequaternion,longitude,latitude,height,andeast,north,verticalvelocitiesofthespaceship.(2) drawthelatitude-versus-longitudetrajectoryofthespaceship,withhorizontallongitudeaxis.(3) drawthecurveoft

48、heheightofthespaceship,withhorizontaltimeaxis.(4) drawthecurvesoftheattitudeanglesofthespaceship,withhorizontaltimeaxis.2.Procedurecode2.Procedurecode.quaternionmultiplycode:quaternionmultiplycode:functionq1=quml(q1,q2);lm=q1(1);p1=q1(2);p2=q1(3);p3=q1(4);q1=lm-p1-p2-p3;p1lm-p3p2;p2p3lm-p1;p3

49、-p2p1lm*q2;end.quaternioninversioncode:quaternioninversioncode:functionqni=qinv(q)q(1)=q(1);q(2)=-q(2);q(3)=-q(3);q(4)=-q(4);qni=q;end. case1DCMalgorithmcase1DCMalgorithm：functionans11cz=cos(80/180*pi)sin(80/180*pi)0;-sin(80/180*pi)cos(80/180*pi)0;001;cx=100;0cos(90/180*pi)sin(90/180*p

50、i);0-sin(90/180*pi)cos(90/180*pi);cy=cos(170/180*pi)0-sin(170/180*pi);010;sin(170/180*pi)0cos(170/180*pi);A=cz*cy*cx*0;0;-9.79end. case1quaternionalgorithmcase1quaternionalgorithmfunctionans12q1=cos(40/180*pi);sin(40/180*pi);0;0;%Thefirstrotationquaternionq2=cos(45/180*pi);0;sin(45/180*pi);0;

51、%Thesecondrotationquaternionq3=cos(85/180*pi);0;0;sin(85/180*pi);%callthequaternionmultiplicationsubfunctionq=quml(r,q3)P1=q(1)q(2)q(3)q(4);-q(2)q(1)q(4)-q(3);-q(3)-q(4)q(1)q(2);-q(4)q(3)-q(2)q(1);P2=q(1)-q(2)-q(3)-q(4);q(2)q(1)q(4)-q(3);q(3)-q(4)q(1)q(2);q(4)q(3)-q(2)q(1);P=P1*P2;gn=P*g;A=gn(2:4)en

52、d%ThethirdrotationDCM%ThefirstrotationDCM%ThesecondrotationDCMg=0;0;0;-9.79;r=quml(q1,q2);. case2SINSquaternionalgorithmcodecase2SINSquaternionalgorithmcodefunctionans2clc;clear;loadmission.matR=6371000;g=9.79;vie=7.292e-5;wie=7.292e-5;lambda=zeros(1,9000);eul=zeros(9000,3);eul(1,:)=-90900;ph

53、i=lambda;h=lambda;VE=lambda;VII=lambda;VI=lambda;fK=lambda;fy=lambda;fz=lambda;avVE=lambda;avVN=lambda;avVT=lambda;lambda(1)=110+57/60;phi(1)=19+37/60;h(1)=100;Q=quatmultiply(cos(-pi/4)00sin(-pi/4),cos(pi/4)sin(pi/4)00);VE(1)=0;VN(1)=0;VT(1)=0;T=1800;dett=0.01;detT=0.2;K=T/detT;k=detT/dett;forN=1:Kf

54、orn=1:kifnorm(GGM(N-1)*k+n,:)=0deltaTheta0=pi*0.01*norm(GGM(N-1)*k+n,:)/(3600*180);DELTATHETA=0.01*pi*.0-GGM(N-1)*k+n,1)-GGM(N-1)*k+n,2)-GGM(N-1)*k+n,3);GGM(N-1)*k+n,1)0GGM(N-1)*k+n,3)-GGM(N-1)*k+n,2);GGM(N-1)*k+n,2)-GGM(N-1)*k+n,3)0GGM(N-1)*k+n,1);GGM(N-1)*k+n,3)GGM(N-1)*k+n,2)-GGM(N-1)*k+n,1)0/(36

55、00*180);Q=(1-deltaTheta0A2/8+deltaTheta0A4/384)*eye(4)+(1/2-deltaTheta0/48)*DELTATHETA)*Q);endendwig=-VN(N)/(R+h(N)VE(N)/(R+h(N)+wie*cosd(phi(N)VE(N)/(R+h(N)*tand(phi(N)+wie*sind(phi(N);detr=wig*detT;detr0=norm(detr);n=detr/detr0;qGrGO=cos(detr0/2)sin(detr0/2)*n;invq=quatconj(qGrGO);Q=quatmultiply(i

56、nvq,Q);eul(N+1,:)=rad2deg(q2eul(Q);invQ=quatinv(Q);fx(N)=AAM(N,1)*1e-7*g;fy(N)=AAM(N,2)*1e-7*g;fz(N)=AAM(N,3)*1e-7*g;fb=fx(N)fy(N)fz(N);fg=quatrotate(invQ,fb);Wie=0wie*cosd(phi(N)wie*sind(phi(N);AGr=fg-cross(wig+Wie),VE(N)VN(N)VT(N)-00g*(RA2/(R+h(N)A2);VE(N+1)=AGr(1)*detT+VE(N);VN(N+1)=AGr(2)*detT+V

57、N(N);VT(N+1)=AGr(3)*detT+VI(N);avVE(N)=(VE(N+1)+VE(N)/2;avVN(N)=(VN(N+1)+VN(N)/2;avVT(N)=(VT(N+1)+VT(N)/2;lambda(N+1)=detT*avVE(N)/R/cosd(phi(N)*180/pi+lambda(N);phi(N+1)=detT*avVN(N)/R*180/pi+phi(N);h(N+1)=detT*avVT(N)+h(N);endattitudequaternion=Qlongitude=lambda(N)latitude=phi(N)height=h(N)eastvel

58、ocity=VE(N)northvelocity=VN(N)verticalvelocity=VI(N)startpoint=R*cosd(phi(1)*sind(lambda(1)-cosd(lambda(1)endpoint=R*cosd(phi(K)*sind(lambda(K)-cosd(lambda(K)R*sind(phi(K);R*sind(phi(1);changeangle=acos(startpoint*endpoint/norm(startpoint)/norm(endpoint);diatanceSE=R*changeanglefigure;plot(lambda,phi)ylabel(latitude/?)xlabel(longitude/?)title(pathoftherocketfigu