Exploitation des observations Obtenues avec un sky Brightness monitor au dôme C  




Annexe 3 : Programmes




Presentation d'un programme IDL puis du programme C que nous avons réalisé.



Programme IDL (sur la donnée du 2 Janvier 2009, 7h30)


Le programme suivant permet de calculer la brillance du ciel autour du soleil en fonction de la distance depuis le centre du soleil.
Un exemple de résultat que génère ce programme est ici.



pro distance
 
dcs = fltarr(320,240)
dcss = dcs
dcs2 = dcs
dcs3 = dcs
moydat= fltarr(70,4)

data = readfits('/home/jarnaud/sbm/usbdisk/dec0827jan0903/sbm6.20090102.0730')

for i = 0,319 do begin & for j=0,239 do begin & dcs(i,j) = sqrt((i-170.)^2+(j-145.)^2) & endfor & endfor
dcss = dcs
for i = 0,168 do begin & for j=0,150 do begin dcs2(i,j) = 2000.  & endfor & endfor
for i = 0,168 do begin
for j=0,150 do begin
dcs3(i,j)=  (float(j)-145.)/dcs(i,j)
endfor & endfor

sunbr = fltarr(4)
for im = 0,3 do begin
sunbr(im) = mean(data(168:178,135:145,im))
dat1 = data(*,*,im)/sunbr(im)

for dd= 9,210,3 do begin
dmi = dd-1
dma=dd+1
dcss(where(dcs ge dmi and dcs le dma)) = 300.
dcss(where(fix(dcs) eq dd)) = 300.
datsec= where(dcss eq 300. and dcs2 eq 2000.)

dcss(datsec) = 1000.
moydat(dd/3-3,im) = mean(dat1(datsec))
;window,1, xs = 320, ys = 240, retain=2
;tvscl, dcss

dcss = fltarr(320,240)
E:
endfor
endfor
window,2, retain=2
; to plot in Rsun units
xgr = findgen(66)
wgra=xgr/7.7
; a very crude correction of scettered light
moydat = moydat -0.15
; Plots of brigthness in 10(-6) Bsun units
plot, wgra, 10*moydat(0:65,0)/moydat(0,0), yr=[0,10], background = 255 , color = 0,  /xst $
, xtitle = 'Distance to the Sun center in units of solar radius', ytitle = 'B in 10(-6) BSun units' $
, title = ' Sky brightness for the 450 nm, 530 nm, 890 nm and 940 nm filters'
oplot, wgra, 10*moydat(0:65,1)/moydat(0,1), color = 0,l=2
oplot, wgra, 10*moydat(0:65,2)/moydat(0,2), color = 0
oplot, wgra, 10*moydat(0:65,3)/moydat(0,3), color = 0,l=2
write_bmp, 'Bsundist200901020730.bmp', tvrd()
stop

end





Haut de page




Programme C


Le programme que nous avons réalisé dans le langage C permet de calculer les coordonnées du soleil en fonction de la position, de la date et et de l'heure. Il permet également de tracer la trajectoire de celui-ci dans le ciel.




#include <stdio.h>
#include <math.h>
#include <stdlib.h>
/*#include "astro.h"*/

//astro.c
//liste des fonctions et routines de calculs utilisés en astro

//------------------------------------------------
//structure necessaire pour les coordonnées equatorial
//------------------------------------------------


#define PI180 0.01745329252
#define N 48

//structure utile pour stocker les données
typedef struct coordonne{
        double xeq;
        double yeq;
        double zeq;
        }coor;

typedef struct coordonne_sol{
        double asc;
        double decl;
        }coorsol;

//---------------------------------------------------
//jour julien
//---------------------------------------------------

double jj(double a,double m,double j,double heure)
{
          double jjj;
          jjj=int(365.25*(a+4716))+int(30.6001*(m+1))+j+heure/24+int(heure)/1440+2-int(a/100)+int(a/400)-1524.5;
          return (jjj);
}

// anomalie moyenne
double anmoy(double a,double m,double j,double heure)
{
       double jjj,M;
       double J2000=2451545.0;
       double M0=357.5291;
       double M1=0.98560028;
       jjj=jj(a,m,j,heure);

       M=M0+M1*(jjj-J2000);
       return M;

}

//anomalie vrai
double anvrai(double M)
{
       double C,nu;
       C=1.9148*sin(M*PI180)+0.0200*sin(2*M*PI180)+0.0003*sin(3*M*PI180);
       nu=C+M;
       return nu;
}

//coordonnée du soleil dependant de l'anomalie moyenne et vrai
coorsol Cs(double nu)
{
       double lambda_s;
       coorsol A;
       double EPS=23.45;
       lambda_s=nu+102.9372+180;
       A.asc=atan(tan(lambda_s*PI180)*cos(EPS*PI180));
       A.decl=asin(sin(lambda_s*PI180)*sin(EPS*PI180));
      /* if (A.asc<0)
          A.asc= -A.asc;*/
       return A;
}

//calcul de l'azimut
double az(double latit,double d,double H){
       double Azimut;
       Azimut=atan(sin(H*PI180)/((sin(latit*PI180)*cos(H*PI180))-(cos(latit*PI180)*tan(d))));
       return Azimut;
}
// coordonnée du soleil dans le référentiel locale
coorsol hh(double latit,double d,double A,double H)
{
       coorsol h;
       h.asc=acos((cos(d)*sin(H/PI180))/(sin(A)));
       h.decl=asin(cos(latit*PI180)*cos(d)*cos(H*PI180)+sin(latit*PI180)*sin(d));
       return h;
}


//fonction principale
int main(){

    double longit,latit,heure;
    double M[N],nu[N];
    double Ma[N];
    double azimut[N],t[N];
    double a,m,j,heure_loc[N],H;
    coor B,S;
    coorsol h[N];
    //coorsol h_vrai;
    int i;
    FILE* pf1;
    FILE* pf2;
    FILE* pf3;
    coorsol A[N];

    pf1=fopen("var_sol","w");
    pf2=fopen("sol_ob","w");
    pf3=fopen("masse_air","w");

    longit=123;
    latit=-75;
    a=2009;
    m=12;
    j=21;


    for (i=0;i<N;i++){
        A[i].asc=0;
        A[i].decl=0;
        H=180; //H a minuit t=12 à longitude =0
        t[i]=0;
        t[i]=jj(a,m,j,i);
        t[i]-=longit/360;//heure a longitude donnée
        heure_loc[i]=(t[i]-(int)(t[i]))*24;
        H+=heure_loc[i]*15;




        M[i]=anmoy(a,m,j,heure_loc[i]);
        nu[i]=anvrai(M[i]);
        A[i]=Cs(nu[i]);
        azimut[i]=az(latit,A[i].decl,H);
        h[i]=hh(latit,A[i].decl,azimut[i],H);
        M[i]=-6378.14*sin(h[i].decl)+sqrt(6378.14*6378.14*sin(h[i].decl)*sin(h[i].decl)+2*6378.14*100+100*100);
        fprintf(pf1,"%lf\t%lf\n",A[i].asc*24/6.28+12,A[i].decl/PI180);

        fprintf(pf2,"%d\t%lf\n",i,h[i].decl/PI180);
        fprintf(pf3,"%lf\t%lf\n",M[i],h[i].decl/PI180);

    }

    fclose(pf1);
    fclose(pf2);




    system("PAUSE");
    return 0;
}






Haut de page