Rambler's Top100Astronet    
  po tekstam   po forumu  vnutri temy
 

args[0]=message
args[1]=DB::DB::Message=HASH(0x3538cc0)
Re: Formuly sfericheskoi trigonometrii tochnye ili priblizhennye ?
2.09.2010 8:29 | S. Yu. Yudin

(povtor)

Zanimayas' voprosom efemeridnoi popravki ET-UT dlya sovremennyh teorii i dlya teorii drevnih astronomov, ya sozdal eshe odnu formu dlya svoei programmy Solsys7 (skrinshot smotrite nizhe), gde u menya vosproizvoditsya animaciya solnechnyh i lunnyh zatmenii i tranzita planet po disku Solnca. Proveryaya vremya nastupleniya zatmenii ili tranzita planet, ya zametil, chto v ekvatorial'nyh koordinatah po animacionnoi kartinke (nizhnyaya kartinka na skrinshote) poluchayutsya horoshie rezul'taty, kak dlya sovremennyh nablyudenii, tak i dlya drevnih, a v gorizontal'nyh koordinatah nablyudaetsya kakaya to strannost'. Bol'shinstvo nablyudenii tozhe horosho vosproizvodilos' i v gorizontal'nyh koordinatah, no inogda nablyudeniya ochen' sil'no rashodilis' s raschetnymi dannymi, kotorye ya poluchal, nablyudaya animacionnuyu kartinku zatmeniya ili tranzita.


http://ser.t-k.ru/Ris/EQUHOR1.jpg

Skrinshot formy 16 programmy Solsys7, gde vosproizvoditsya tranzit Venery 8.6.2004 (na nizhnei kartinke animaciya v gorizontal'nyh koordinatah), gde tochkami 07 oboznacheny vosem' tochek na diske Solnca, koordinaty kotoryh opredeleny v ekvatorial'noi sisteme koordinat i zatem pereschitany v gorizontal'nye koordinaty. Kruzhok na diske Solnca eto Venera.


V centre kartinki ya risoval vidimyi diametr Solnca (dlya solnechnyh zatmenii i tranzitov), a zatem po raznosti koordinat (ekvatorial'nyh ili gorizontal'nyh) mezhdu centrom diska Solnca i planety ili Luny nahodil polozhenie na kartinke centra planety ili Luny i risoval ih vidimyi diametr. Pri etom, pervonachal'no ya poluchal dlya sovremennyh teorii N'yukoma ili JPL (Laboratoriya reaktivnogo dvizheniya podrazdelenie NASA) geliocentricheskie eklipticheskie koordinaty planet (ili v epohe daty ili v standartnoi epohe J2000) dlya planet i geocentricheskie eklipticheskie koordinaty Luny po teorii Brauna (v epohe daty), a zatem privodil ih k vidimym koordinatam, t.e. k nablyudaemym s konkretnoi tochki Zemli. Analogichnye preobrazovaniya ya delal i s raschetnymi dannymi, poluchennymi po tablicam drevnih astronomov (Ptolemei, Al' Horezmi, Kepler i t.d.). Na skrinshote v okoshkah sprava otrazhena posledovatel'nost' preobrazovanii. Cepochka preobrazovanii poluchaetsya dovol'no taki dlinnoi po etomu v pervuyu ochered' ya, konechno zhe, zapodozril kakuyu to oshibku v kode programmy, no proverka (pryamoe preobrazovanie koordinat i obratnoe) pokazala, chto oshibki net.


Togda ya zanyalsya proverkoi standartnyh funkcii preobrazovaniya, kotorye ya vzyal iz knigi O. Montenbruk i T. Pfleger Astronomiya s personal'nym komp'yuterom i perevel ih s yazyka programmirovaniya Pascal na yazyk Visual Basic 6.0, a v pervuyu ochered' ya zanyalsya funkciei perevoda ekvatorial'nyh koordinat v gorizontal'nye i neobhodimoi dlya etogo funkciei rascheta Grinvichskogo zvezdnogo vremeni GMST (v kode funkcii posle apostrofa ` idet kommentarii k kodu).


Private Sub EQU_HOR(ByVal j As Integer)
TAU = ObsLMST - GeoL(j): If TAU < 0 Then TAU = TAU + 360 'chasovoi ugol v gradusah
Call EQUHOR(j) 'funkciya perevoda koordinat iz knigi O. Montenbruk i T. Pfleger
Call AlfaBetta(15, j) 'nahodim po dekartovym koordinatam polyarnye koordinaty
Azimut(j) = AlfaGR(j): If Azimut(j) > 180 Then Azimut(j) = Azimut(j) - 360
Visota(j) = Betta1
End Sub


Public Sub EQUHOR(i As Integer) 'conversion of equatorial into horizontal coordinates
'analog PROCEDURE EQUHOR (DEC,TAU,PHI: REAL; VAR H,AZ: REAL);
'(* DEC = GeoB(I) : declination (-90 deg .. +90 deg) *)
'(* TAU : hour angle (0 deg .. 360 deg) *)
'(* PHI = ObsB : geographical latitude (in deg) *)
Dim CS_PHI#, SN_PHI#, CS_DEC#, SN_DEC#, CS_TAU#
CS_PHI = Cos(ObsB / kUgol): SN_PHI = Sin(ObsB / kUgol)
CS_DEC = Cos(GeoB(i) / kUgol): SN_DEC = Sin(GeoB(i) / kUgol)
CS_TAU = Cos(TAU / kUgol)
X(i) = CS_DEC * SN_PHI * CS_TAU - SN_DEC * CS_PHI
Y(i) = CS_DEC * Sin(TAU / kUgol)
Z(i) = CS_DEC * CS_PHI * CS_TAU + SN_DEC * SN_PHI
'POLAR (X,Y,Z, DUMMY,H,AZ)
End Sub


Public Function LMST(MJD As Double) As Double ' rasschityvaetsya GMST v chasah
'FUNCTION LMST(MJD,LAMBDA:REAL):REAL; local mean sidereal time
'MJD = JD - 2400000.5
T2 = (MJD - Fix(MJD)) * 24: T1 = (Fix(MJD) - 51544.5) / 36525
LMST = 6.697374558 + 1.0027379093 * T2 + (8640184.812866 + (0.093104 - 0.0000062 * T1) * T1) * T1 / 3600# End Function

K szhaleniyu, proverka pokazala, chto funkcii dayut te zhe rezul'taty, chto i primery opisannye v literature, naprimer, v knige Jean Meeus Astronomical algorithmus. V knige dayutsya raschetnye topocentricheskie ekvatorial'nye koordinaty Venery na 19:21:00 UT dlya observatorii USNO (ObsL= -77,0656 ObsB= 38,9214) - dolgota 347,3193 gradusov i shirota -6,72 gradusa pri GMST=08:34:56,853. A u menya programma daet raschetnye ekvatorial'nye koordinaty dolgota 347,3168 i shirota -6,7216 pri GMST=08:34:57 (sekundy okrugleny i GMST poschitano bez ucheta nutacionnoi popravki). Pri etih ishodnyh dannyh v knige poluchayutsya sleduyushie gorizontal'nye koordinaty azimut 68,0337 i vysota 15,1249, a u menya programma daet azimut 68,035 i vysota 15,1211. Rashozhdeniya poluchilis' v pyatom znake, chto vpolne dopustimo pri nebol'shih razlichiyah v ishodnyh ekvatorial'nyh koordinatah i pri otsutstvie u menya nutacionnoi popravki pri raschete GMST.


Voobshe-to provedennaya mnoyu proverka byla prosto shagom otchayaniya, t.k. vse eti funkcii ya i tak mnogokratno testiroval na drugih formah programmy Solsys i nikakih zamechanii oni u menya ne vyzyvali. No tak kak poluchennye mnoyu rezul'taty v ekvatorial'nyh i gorizontal'nyh koordinatah sil'no otlichalis' (chego ne dolzhno byt' v principe) ya reshil vospol'zovat'sya eshe i formulami perevoda ekvatorial'nyh koordinat v gorizontal'nye i v polyarnyh koordinatah (privedennaya vyshe procedura EQUHOR proizvodit perevod koordinat v dekartovoi sisteme koordinat, a potom vychislyayutsya polyarnye koordinaty). Ya vzyal standartnye funkcii perevoda (privodyatsya v lyubom uchebnike) i napisal proceduru EQU_HOR2. Pri etom, t.k. v yazyke programmirovaniya Visual Basic 6.0 net mnogih trigonometricheskih funkcii (kak i vo mnogih drugih yazykah), naprimer, nuzhnyh mne arcsin i arccos, a est' tol'ko arctan (v kode Atn), to mne prishlos' napisat' dlya nih kod samomu (ispol'zuya izvestnye trigonometricheskie preobrazovaniya).


Private Sub EQU_HOR2(ByVal j As Integer)
TAU = ObsLMST - GeoL(j): If TAU < 0 Then TAU = TAU + 360 ' TAU chasovoi ugol v gradusah
'H=sin(Visota(j)) vremennaya peremennaya
H = Sin(GeoB(j) / kUgol) * Sin(ObsB / kUgol) + Cos(GeoB(j) / kUgol) * Cos(ObsB / kUgol) * Cos(TAU / kUgol)
Visota(j) = Atn(H / (1 - H * H) ^ 0.5) * kUgol ' funkciya arcsin -1...+1
H = Sin(GeoB(j) / kUgol) - Sin(ObsB / kUgol) * Sin(Visota(j) / kUgol)
H = H / Cos(ObsB / kUgol) / Cos(Visota(j) / kUgol) ' H=cos(Azimut(j)) vremennaya peremennaya
Azimut(j) = Atn((1 - H * H) ^ 0.5 / H) * kUgol ' funkciya arccos 0...+1
AZ = Azimut(j)
If AZ >= 0 And TAU <= 180 Then Azimut(j) = 180 - AZ
If AZ >= 0 And TAU > 180 Then Azimut(j) = AZ - 180
If AZ < 0 And TAU <= 90 Then Azimut(j) = -AZ
If AZ < 0 And TAU > 90 Then Azimut(j) = AZ
End Sub


Rezul'taty poluchilis' te zhe samye, chto i pri preobrazovaniyah v dekartovoi sisteme koordinat, i ya byl prosto v tupike ne znaya, chto mozhno eshe proverit'. Togda ya reshil prosto tupo posmotret', a kak perevodyatsya iz ekvatorial'nyh koordinat v gorizontal'nye ne tol'ko koordinaty centra Solnca, no eshe i vos'mi tochek raspolozhennyh na diske Solnca. I tut okazalos', chto v kakih to polozheniyah Solnca na gorizonte tochki raspolozhennye po okruzhnosti Solnca v ekvatorial'nyh koordinatah perevodyatsya v gorizontal'nye pravil'no i my vidim tu zhe okruzhnost', a v kakih to polozheniyah my vidim v gorizontal'nyh koordinatah ellips. Naprimer, na skrinshote, pokazan vid Solnca pri azimute 0 gradusov (napravlenie na Yug), a pri azimutah -90, +90 i +/-180 gradusov priveden na risunke nizhe.


http://ser.t-k.ru/Ris/EQUHOR2.jpg


Kak vidim pri azimute +/-180 gradusov koordinaty tochek perevodyatsya tak, chto sovpadayut s diskom Solnca narisovannom otnositel'no centra Solnca, koordinaty kotorogo tozhe pereschitany v gorizontal'nye koordinaty, a pri azimutah -90 i +90 gradusov u nas koordinaty vos'mi tochek neznachitel'no otstoyat ot okruzhnosti, t.e. perevodyatsya pochti bez oshibki. No v raznoe vremya goda maksimal'naya oshibka pri perevode koordinat poluchaetsya pri raznyh azimutah i, naprimer, v dekabre ona budet maksimal'noi uzhe ne pri 0 gradusov kak letom, a pri azimute +/-180, a v marte i sentyabre i pri 0 gradusov i pri +/-180 gradusov. Vot tol'ko pri etom velichina maksimal'noi oshibki po vidimomu radiusu Solnca v marte i sentyabre budet tol'ko 0,1 gradusa, a v iyune i dekabre uzhe pochti 0,4 gradusa, chto ochen' mnogo i sopostavimo s vidimym diametrom Solnca (0,5 gradusa). I esli my budem obrabatyvat' dannye nablyudenii za Solncem, kotorye vypolneny ne po centru diska, a po kakomu-to ego krayu, to my poluchim ochen' ne tochnye dannye.


Takim obrazom, ya ponyal otkuda u menya poluchayutsya raznye dannye v raznyh sistemah koordinat dlya zatmenii i tranzitov planet i voznikaet vopros. A mozhet byt' formuly sfericheskoi trigonometrii yavlyayutsya ne tochnymi, a priblizhennymi i pri opredelennyh uglah dayut znachitel'nuyu pogreshnost', chto ya i poluchil. Vot otvet na etot vopros ya by i hotel poluchit'. Ili mozhet byt' kto-to znaet bolee tochnye formuly perevoda ekvatorial'nyh koordinat v gorizontal'nye. A esli kto-to schitaet, chto u menya v programme imeetsya kakoi to glyuk, to u menya k nemu budet pros'ba vzyat' ekvatorial'nye koordinaty 8-i tochek po okruzhnosti Solnca 22 iyunya v 12 chasov mestnogo vremeni dlya punkta nablyudeniya nahodyashegosya v severnom polusharii, perevesti ih svoimi programmami ili obychnym raschetom v gorizontal'nye i nanesti ih na risunok.

Kstati, zamechennyi mnoyu effekt ellipsnosti Solnca v gorizontal'nyh koordinatah mozhno nablyudat' i v kul'tovoi programme StarCalc (u menya versiya 5.7) esli otslezhivat' dvizhenie Solnca pri ochen' bol'shom masshtabe (20000%40000%), kak eto pokazano na skrinshote nizhe. Pri etom v nastroikah dlya koordinat luchshe vybrat' dlya izobrazheniya neba ne parallel'nuyu, a konicheskuyu setku koordinat, t.k. pri etom na animacionnoi kartinke Solnce budet risovat'sya v vide okruzhnosti, a ego ellipsnost' mozhno budet opredelit' po koordinatnoi setke (pri parallel'noi proekcii budet i neravnomernost' masshtaba koordinatnoi setki i ellipsnost' Solnca, chto smazhet etot effekt). Na skrinshote dazhe vizual'no vidno, chto po koordinatnoi setke razmer po azimutu primerno v dva raza bol'she chem po vysote, a konkretno po vysote razmer budet 67 gradusov 17 minut 66 gradusov 45 minut = 32 minuty (real'nyi vidimyi diametr Solnca), a po azimutu 0 gradusov 40 minut (-0 gradusov 40 minut) = 80 minut. V rezul'tate poluchaem, chto oshibka po vidimomu radiusu Solnca poluchaetsya (80-32)/2=24 minuty, t.e. te zhe samye 0,4 gradusa, chto ya poluchil po animacionnoi kartinke programmy Solsys7.


http://ser.t-k.ru/Ris/EQUHOR3.jpg


Takim obrazom, my eshe raz ubezhdaemsya, chto pri perevode ekvatorial'nyh koordinat neskol'kih tochek raspolozhennyh po okruzhnosti Solnca v gorizontal'nye my poluchaem ne okruzhnost', a ellips, chego ne dolzhno byt' v principe, t.k. pri prostom povorote sistemy koordinat (masshtaby po osyam ostalis' neizmennymi) takoi effekt ne dolzhen nablyudat'sya. No, esli u Vas pri takom perevode koordinat poluchitsya okruzhnost', to Vy obyazatel'no soobshite mne kakimi formulami Vy pol'zovalis'. A, esli kakaya to oshibka imeetsya v moih rassuzhdeniyah, to ukazhite gde u menya oshibka v rassuzhdeniyah, po tomu, chto ya, naprimer, ne ponimayu pochemu programma StarCalc pri ispol'zovanii funkcii opredeleniya uglovogo rasstoyaniya mezhdu dvumya tochkami daet po azimutu rezul'tat 32 minuty, hotya po koordinatnoi setke mezhdu etimi zhe dvumya tochkami my poluchili 80 minut. Hotya, byt' mozhet, takoi problemy s formulami sfericheskoi trigonometrii i ne sushestvuet i prosto ne hvataet tochnosti vychislenii komp'yutera. Ya, naprimer, ispol'zuyu v programme 16 razryadov, a v programme StarCalc, kak ya dogadyvayus', ispol'zuetsya 32 razryada, no po grubomu prikidu dlya ustraneniya oshibki v desyatyh dolyah gradusa hvatilo by i 8 razryadov (ili ya ne prav).

S nailuchshimi pozhelaniyami Sergei Yudin.




Forumy >> Astronomiya i Internet
Spisok  /  Derevo
Zagolovki  /  Annotacii  /  Tekst

Astronet | Nauchnaya set' | GAISh MGU | Poisk po MGU | O proekte | Avtoram

Kommentarii, voprosy? Pishite: info@astronet.ru ili syuda

Rambler's Top100 Yandeks citirovaniya