Des grumeaux dans la trigo (Pascal vs VBA)

Version imprimable

Bonjour :)
une nouvelle tracasserie à vous proposer
je réalise la motorisation de mon télescope
moteurs pas à pas + arduino +pc
pour fonctionner, le soft a besoin de connaitre l'altitude et l'azimut de l'objet pour en assurer le suivi
on m'a passé un bout de code en VB que j'ai transcrit en Pascal mais j'ai dû me tromper car les fonctions ne renvoient pas les valeurs correctes
les fonctions VB

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
Function Angle_horaire(AD As Date, Longitude As Single, Date_choisie As Date, Heure_choisie As Date, Zone As Integer) As Single
 
Dim Dayz, Monthz, Yearz, Hourz, Minutez, Secondz As Integer
Dim AD_in_degre, A, B, C, D, JJ, T, H1, HSH, HS, AngleH, AngleT, PI As Single
 
PI = 3.14159
Dayz = Day(Date_choisie)
Monthz = Month(Date_choisie)
Yearz = Year(Date_choisie)
Hourz = Hour(Heure_choisie)
Minutez = Minute(Heure_choisie)
Secondz = Second(Heure_choisie)
 
AD_in_degre = 15 * (Hour(AD) + Minute(AD) / 60 + Second(AD) / 3600)
 
If (Monthz < 3) Then
Monthz = Monthz + 12
Yearz = Yearz - 1
End If
 
A = Int(Yearz / 100)
B = 2 - A + Int(A / 4)
C = Int(365.25 * Yearz)
D = Int(30.6001 * (Monthz + 1))
JJ = B + C + D + Dayz + 1720994.5
T = (JJ - 2451545) / 36525
H1 = 24110.54841 + (8640184.812866 * T) + (0.093104 * (T * T)) - (0.0000062 * (T * T * T))
HSH = H1 / 3600
HS = ((HSH / 24) - Int(HSH / 24)) * 24
AngleH = (2 * PI * HS / (23 + 56 / 60 + 4 / 3600)) * 180 / PI
AngleT = ((Hourz - 12 + Minutez / 60 - Zone) * 2 * PI / (23 + 56 / 60 + 4 / 3600)) * 180 / PI
Angle_horaire = AngleH + AngleT - AD_in_degre + Longitude
 
End Function
 
Function Convert_temps(ByVal T As Single) As String
Dim H, M, S As Integer
Dim Temp, Temp2 As Single
 
H = Int(T)
If H >= 24 Then H = H - 24
Temp = T - H
Temp2 = Temp * 60
M = Int(Temp2)
If M >= 60 Then
M = 0
H = H + 1
End If
 
Temp = Temp - M / 60
Temp2 = Temp * 3600
S = Temp2
If S >= 60 Then
S = 0
M = M + 1
End If
 
Convert_temps = H & ":" & M & ":" & S
 
End Function
 
Function Arsin(ByVal X As Single) As Single
If (Abs(X) >= 1) Then Arsin = 0 Else Arsin = Atn(X / Sqr(-X * X + 1))
End Function
 
Function Arcos(ByVal X As Single) As Single
If (Abs(X) >= 1) Then Arcos = 0 Else Arcos = Atn(-X / Sqr(-X * X + 1)) + 2 * Atn(1)
End Function
 
Function Hauteur(Dec As Single, Latitude As Single, H As Single) As Single
 
Dim Sin_hauteur, PI As Single
PI = 3.14159
Sin_hauteur = Sin(Dec * PI / 180) * Sin(Latitude * PI / 180) - Cos(Dec * PI / 180) * Cos(Latitude * PI / 180) * Cos(H * PI / 180)
Hauteur = Arsin(Sin_hauteur) * 180 / PI
End Function
 
Function Azimuth(Dec As Single, Lat As Single, H As Single, Haut As Single) As Single
 
Dim Cosazimuth, Sinazimuth, test, Az As Single
PI = 3.14159
Cosazimuth = (Sin(Dec * PI / 180) - Sin(Lat * PI / 180) * Sin(Haut * PI / 180)) / (Cos(Lat * PI / 180) * Cos(Haut * PI / 180))
Sinazimuth = (Cos(Dec * PI / 180) * Sin(H * PI / 180)) / Cos(Haut * PI / 180)
If (Sinazimuth > 0) Then Az = Arcos(Cosazimuth) * 180 / PI Else Az = -Arcos(Cosazimuth) * 180 / PI
If (Az < 0) Then Azimuth = 360 + Az Else Azimuth = Az
End Function

le Pascal

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
Function TForm1.Angle_horaire(AD : tdate; Longitude : Single; Date_choisie :
tdate; Heure_choisie : tdate; Zone : Integer) : Single  ;
var
  Date1   : TdateTime;
  Annee, Mois, Jour, Heure, Minute, Seconde, milliSec,xannee,xmois,xjour,xheure,xminute,xseconde,xmillisec : Word;
  NbJour  : Integer;
  JourJul : Double;
  AD_in_degre, A, B, C, D, JJ, T, H1, HSH, HS, AngleH, AngleT, PI : Single  ;
begin
  PI := 3.14159 ;
  Date1 := now;
  decodeDateTime(now, Annee, Mois, Jour, Heure, Minute, Seconde, milliSec);
  decodeDateTime(ad,xannee,xmois,xjour, xheure, xminute, xseconde,xmillisec);
  JourJul := DateTimeToJulianDate(Date1);// + NbJour;
AD_in_degre := 15 * (xheure + (xminute / 60) + (xseconde / 3600));
  If (mois < 3) Then begin
  mois := mois + 12  ;
  annee := annee - 1 ;
  End ;
 A := Int(annee / 100);
B := 2 - A + Int(A / 4) ;
C := Int(365.25 * annee) ;
D := Int(30.6001 * (mois + 1)) ;
JJ := B + C + D + jour + 1720994.5 ;
T := (JJ - 2451545) / 36525     ;
H1 := 24110.54841 + (8640184.812866 * T) + (0.093104 * (T * T)) -
(0.0000062 * (T * T * T))  ;
HSH := H1 / 3600 ;
HS := ((HSH / 24) - Int(HSH / 24)) * 24  ;                     ;
AngleH := (2 * PI * HS / (23 + 56 / 60 + 4 / 3600)) * 180 / PI  ;
AngleT := ((heure - 12 + (Minute / 60) - Zone) * 2 * PI / (23 + 56 / 60 + 4 /
3600)) * 180 / PI   ;
Angle_horaire := AngleH + AngleT - AD_in_degre + Longitude   ;
end;
function tform1.Arsin (X : Single) : Single   ;
begin
If (Abs(X) >= 1) Then Arsin := 0 Else Arsin := arctan(X / sqr(-X * X + 1))  ;
End ;
Function TForm1.Arcos( X : Single) : Single;
begin
If (Abs(X) >= 1) Then Arcos := 0 Else Arcos := arctan(-X / Sqr(-X * X + 1)) + 2 *
arctan(1);
End ;
Function TForm1.Hauteur(Dec : Single; Latitude : Single; H : Single) : Single  ;
var
Sin_hauteur, PI : Single ;
begin
PI := 3.14159 ;
Sin_hauteur := Sin(Dec * PI / 180) * Sin(Latitude * PI / 180) - Cos(Dec * PI /
180) * Cos(Latitude * PI / 180) * Cos(H * PI / 180) ;
Hauteur := Arsin(Sin_hauteur) * 180 / PI;
End ;
Function TForm1.Azimuth(Dec: Single; Lat: Single; H: Single; Haut: Single):Single  ;
var
Cosazimuth, Sinazimuth, test, Az, PI : Single ;
begin
PI := 3.14159   ;
Cosazimuth := (Sin(Dec * PI / 180) - Sin(Lat * PI / 180) * Sin(Haut * PI /
180)) / (Cos(Lat * PI / 180) * Cos(Haut * PI / 180));
Sinazimuth := (Cos(Dec * PI / 180) * Sin(H * PI / 180)) / Cos(Haut * PI / 180);
If (Sinazimuth > 0) Then Az := Arcos(Cosazimuth) * 180 / PI Else Az := -Arcos(Cosazimuth) * 180 / PI ;
If (Az < 0) Then Azimuth := 360 + Az Else Azimuth := Az;
End ;

la fonction angle horaire fonctionne bien
mais "hauteur" et "azimut" non (en VB tout fonctionne)
ne connaissant rien en trigo, je me demande si les fonctions pareillement nommées en VB/Pascal donnent les mêmes résultats
j'ai pu vérifier "sin" et "cos" qui semblent bonnes
qu'en penses-vous?
merci de m'avoir lu 8-)

1 pièce(s) jointe(s)

Bonjour,

Déjà, dans un premier temps, je te conseille d'utiliser l'unité math afin de récupérer des fonctions prédéfinies et d'éviter des erreurs d'arrondis sur Pi :mouarf:. De même, j'utiliserais volontiers dans le cas d'astronomie des Double au lieu des Single (mais c'est une affaire de goût).

Avec math, on peut écrire des choses comme :
Code:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 implementation uses math; {$R *.lfm} { TForm1 } procedure TForm1.Button1Click(Sender: TObject); var R: Double; begin R := Pi; R := ArcSin(0.5); R := ArcCos(0.10); end;
Pi, Arccos et Arcsin sont par conséquent déjà définis. La seule différence avec le code fourni (en dehors de Arcsin et ArcCos au lieu de Arsin et Arcos) est qu'il déclenche une exception si les valeurs données en paramètre sont supérieures à 1 (dans le code fourni, les fonctions renvoient 0 :weird:).

Au lieu de la multitude de PI / 180, utilise la fonction degtorad (degrés en radians) qui est plus explicite. De même, 180 / PI est rendu par la fonction radtodeg (radians en degrés).

Pour le résultat des fonctions, je te conseille d'utiliser Result au lieu du nom de la fonction, ce qui permettra de produire du code compatible avec Delphi, par exemple.

De même, pourquoi ne pas regrouper les paramètres de même nature de manière à rendre le code plus lisible ?

Exemple :
Code:

1 2 3 4 5 6 7 function TForm1.Hauteur(const Dec, Latitude, H : Single) : Single; var Sin_hauteur: Single; begin Sin_hauteur := Sin(degtorad(Dec)) * Sin(degtorad(Latitude)) - Cos(degtorad(Dec)) * Cos(degtorad(Latitude)) * Cos(degtorad(H)) ; Result := radtodeg(Arcsin(Sin_hauteur)); end ;
La fonction paraît bien plus lisible, non ? C'est ainsi que je me demande si elle est juste, y compris en VBA : elle est régulièrement donnée avec une addition et non une soustraction... :weird:

Pièce jointe 214392

Je regarde la suite bientôt,

Gilles

30/06/2016, 17h30
droggo

Bonjour,

D'autre part, il est utile de définir des constantes pour les valeurs fixes fréquemment utilisées.

Dans ton cas, il y a par exemple PI/180. :)

Et oui, il vaut mieux utiliser des double.
30/06/2016, 17h51
nomdutilisateur

Déjà quelques pistes :)
je vais reprendre ça ce soir
merci :)
30/06/2016, 19h46
nomdutilisateur

déjà les remèdes conseillés (surtout les trucs issus de l'unité math) portent leurs fruits
la fonction "hauteur" marche
il ne me reste plus qu'à vérifier l'azimut
il faut que je prenne des cas extrêmes pour vérifier la solidité du code (longitudes et latitudes, objets de part et d'autres du méridien et de l'horizon)
restera à tester "pour de vrai" mais là, il faudra du temps :)
02/07/2016, 09h10
nomdutilisateur

la suite....:)
l'utilisation d'un tdate pour transmettre l’ascension droite de l'objet visé amenait des erreurs lors de l'encodage et décodage, je m'en suis apperçu en simulant un suivi sur la Polaire, les courbes variaient au lieu de rester identiques
en utilisant un record fait pour, ça va mieux
il me reste à calculer les variations de vitesse des moteurs selon la hauteur et l'azimut de l'objet
c'est pas évident
un ami me conseille d'utiliser une dérivée:
V à t0 = ( - position à t0 + position à t0+30 secondes ) / 30
mes préoccupations s'éloignent de l'utilisation de Lazarus mais, sait-on jamais, vous avez peut-être la formule magique? 8-)
02/07/2016, 15h13
Archimède

pour calculer la vitesse instantanée d'un point M mobile, les logiciels calculent la dérivée du vecteur position par rapport au temps en faisant :
pour l'instant ti
v(Mi)= Mi-1Mi+1/delta(t)= OMi+1-OMi-1/delta(t)=delta(OM)/delta(t) avec delta(t)=ti+1-ti-1 delta(t) tendant vers 0.

concrètement, tu fais la position d'après moins la position d'avant que tu divises par l'intervalle de temps correspondant.

Si tu as par exemple une acquisition toutes les 100 ms, tu divises par 200 ms (soit 100+100 deux intervalles).
les logiciels comme regressi ou le tableur de latispro procèdent je pense de cette manière pour afficher une dérivée.
03/07/2016, 10h55
nomdutilisateur

Citation:

Envoyé par Archimède

Si tu as par exemple une acquisition toutes les 100 ms, tu divises par 200 ms (soit 100+100 deux intervalles).

Merci de ta proposition, elle confirme mes propres conclusions
en fait, il ne s'agit pas d’acquisitions car avec une focale de 2m, la rotation de champs serait rédhibitoire
c'est pour faire essentiellement du visuel mais cela revient au même avec une fréquence de corrections bien lissée pour le plaisir :)
il me reste à finir et tester tout ça, l'été ne sera pas de trop 8-)

NB: je suis abonné à la discut, mon profil contient un mail valide, rien dans les spams et pourtant rien dans la BAL :o
03/07/2016, 11h46
gvasseur58

Citation:

Envoyé par nomdutilisateur

NB: je suis abonné à la discut, mon profil contient un mail valide, rien dans les spams et pourtant rien dans la BAL :o

As-tu bien choisi le bon rythme et/ou la bonne destination des notifications ? Pour cela, vérifie la valeur choisie dans la liste déroulante en face de la case (cochée !) de l'abonnement.

Cordialement,

Gilles

PS : les notifications fonctionnent en général très bien. Si le problème persiste, me contacter en MP.
03/07/2016, 12h25
nomdutilisateur

Citation:

Envoyé par gvasseur58

As-tu bien choisi le bon rythme et/ou la bonne destination des notifications ? Pour cela, vérifie la valeur choisie dans la liste déroulante en face de la case (cochée !) de l'abonnement.

:oops: c'était ça! merci:D