/* Simulation de donnees mission, Eric Aubourg & Jacques Delabrouille, mai 1999 */ //#include //#include #include #include #include "manip.h" #include "choix_acquisition.h" #include "archeops.h" #include "arcunit.h" #include "choix_param.h" #include "structure.h" #include "tm.h" #include "simulmission.h" #include "simulstate.h" #include "compress.h" #define M_PI 3.14159265358979 double GauRnd(double am, double s); #define drand01() ( (double)rand() / RAND_MAX ) static int numblock = 0; static int stateChanged = 0; static long tickStart = -1; double overSpeed = 1; double whiteNoise[6]; double bolfreq[6]; double glitchFreq; // Hz double glitchMaxAmpl; // W double boloTimeCst; // s int sst2Bars; param_bolo simParam; reglage_bolo simReglage; // Fichier avec donnees mission static short fmisdat = 0; struct boloobs { double power; }; struct posbal { double lon; // degres double lat; // degres double azimut; // degres }; struct misdat { double mjd; // JD - 2450000 --- 9 octobre 1995 midi UTC struct posbal pos; // Position du ballon struct boloobs bolo[6]; // Puissance sur chaque bolo }; struct localgeom { double ra; // Centre du pointage double dec; double dra; // vecteur de deplacement, normalise a une minute d'arc double ddec; double ora; // vecteur orthogonal au deplacement, vers le zenith double odec; // normalise a une minute d'arc }; static int szmisdatbuf; static struct misdat* misdatbuf[2]; static int curBuf = 0; static int iBuf = 0; static int needBuf = -1; static int iRdsample; // celui de lastDat static int iGensample; // 0..oversample-1 static int oversample; // surechantillonnage necessaire... // Donnees GSC pour simulation senseur stellaire struct star { float ra; float dec; short mag; // mag * 100 }; static struct star* stars=0; // trie par dec croissant... long nstars = 0; // Donnees correspondant a un scan long databolo[nb_max_bolo][nb_per_block*2]; int datadiod[nb_per_block*2][48]; double latitude, longitude, mjd; static double per_echant(void); static void restartTime(void); #define lastData misdatbuf[curBuf][iBuf] #define nextData misdatbuf[curBuf][iBuf+1] // Valeur precedente du bolometre pour gestion inertie (glitchs...) en muV static double boloLast[nb_max_bolo]; // Diverses fonctions utiles... void calcLocalGeom(struct misdat* md, struct localgeom* geom); static void restartTime(void) { tickStart = TickCount() - numblock * per_echant() * (double)nb_per_block*2.*60.; } static volatile int hasDebugInfo=0; static volatile char debugInfo[80]; static long readblocks = 0; extern param_bolo parametr; // dans archeops.c extern reglage_bolo reglage_standard[8]; // liste bolo dans le programme principal static double laststarmag[46][10]; static double laststartime[46][10]; void InitSimulMission(void) { long n,nn; long nSamples; for (n=0; n<6; n++) { bolfreq[n] = 100e9; whiteNoise[n] = 2.e-17; boloLast[n] = 0; } for (n=0; n<46; n++) for (nn=0; nn<10; nn++) { laststarmag[n][nn] = laststartime[n][nn] = 0; } glitchFreq = 0.3; glitchMaxAmpl = 1.e-11; boloTimeCst = 1.e-2; // Datacards ReadSimulDC(); // Donnees mission, initialisation du buffer if (FSOpen("\pmission.dat", 0, &fmisdat)) { printf("***ERREUR*** ouverture mission.dat\n"); return; } szmisdatbuf = 512; misdatbuf[0] = (struct misdat*) malloc(szmisdatbuf * sizeof(struct misdat)); misdatbuf[1] = (struct misdat*) malloc(szmisdatbuf * sizeof(struct misdat)); curBuf = 0; iBuf=0; needBuf = -1; iRdsample = 0; iGensample = 0; n = sizeof(long); FSRead(fmisdat, &n, (char*) &nSamples); printf("Fichier mission.dat : %d echantillons\n", nSamples); n = szmisdatbuf * sizeof(struct misdat); FSRead(fmisdat, &n, (char*) misdatbuf[0]); n = (szmisdatbuf-1) * sizeof(struct misdat); memcpy(misdatbuf[1], &misdatbuf[0][szmisdatbuf-1], sizeof(struct misdat)); FSRead(fmisdat, &n, (char*) misdatbuf[1]); readblocks = 2*szmisdatbuf; // Nb d'echantillons dans un cercle, nominalement : 1 cercle = 30 secondes, // ech a 180 Hz, et en realite 125 ech par cercle dans le fichier... oversample = (int) ((30 * 180.) / 125.); printf("oversample = %d\n",oversample); // Catalogue GSC pour senseur stellaire { short ref; long n; //if (FSOpen("\phyp6HADY.dat",0,&ref)) { // printf("***ERREUR*** ouverture hyp6HADY.dat, on va utiliser gsc7\n"); if (FSOpen("\pgsc7.dat",0,&ref)) { printf("***ERREUR*** ouverture gsc7.dat\n"); return; } //} n = sizeof(long); FSRead(ref, &n, (char*) &nstars); n = nstars*sizeof(struct star); stars = (struct star*) malloc(n); FSRead(ref, &n, (char*) stars); FSClose(ref); printf("GSC m<7 : read %d stars\n", nstars); } hasDebugInfo=0; simReglage = reglage_standard[0]; simParam = parametr; } void SimulMissionReadLoop(void) { long n; short rc; if (hasDebugInfo) { printf("*** Debug **%d**\n",hasDebugInfo); hasDebugInfo=0; } if (needBuf < 0) return; if (fmisdat == 0) return; readblocks += szmisdatbuf; printf("Reading simulation data block %d\n",readblocks); memcpy(misdatbuf[needBuf], &misdatbuf[1-needBuf][szmisdatbuf-1], sizeof(struct misdat)); n = (szmisdatbuf-1) * sizeof(struct misdat); if ((rc=FSRead(fmisdat, &n, (char*) misdatbuf[needBuf]))!=0) { printf("Erreur lecture %d\n", rc); FSClose(fmisdat); fmisdat = 0; } needBuf = -1; } float isStarInArray(float ra, float dec, float alow, float dlow, float a1, float d1, float a2, float d2) // -1 = no, 0-1 : yes { float x,y; ra = ra-alow; if (ra>12) ra -=24; if (ra<-12) ra += 24; x = (15*cos(dec * 3.1415926/180)); ra *= x; a1 *= x; a2 *= x; x = ((ra)*a1 + (dec-dlow)*d1)/(a1*a1 + d1*d1); y = ((ra)*a2 + (dec-dlow)*d2)/(a2*a2 + d2*d2); if (x<0 || x>1 || y<0 || y>1) return -1; return y; } // Constantes de temps du senseur stellaire static const double t1 = 17e-3; static const double t2 = 150e-3; static const double t3 = 0; static const double dt = 30e-3; static double sig1(double t); static double sig2(double t); static double sig1(double t) { return (t2*(t1-t3)*exp(-t/t1) - t1*(t2-t3)*exp(-t/t2))/(t1-t2)/t1; } static double sig2(double t) { double x=0; if (t>=0) x -= sig1(t); if (t>=dt) x += sig1(t-dt); return x; } // Balayage dans la direction dra, ddec, pointage central sur ra, dec void SSTSignal(float ra, float dec, float dra, float ddec, float ora, float odec, int* diodes, double secondes) { // Valeurs precedentes pour gestion des constantes de temps //float ora, odec; // orthogonal to dra, ddec, directed to zenith float dmin, dmax; // min and max value of dec in which to look for stars int imin, imax; // corresponding indices in gsc catalog float dlow, alow; // "bottom left" corner of first diode array float d1, a1, d2, a2; // sides of diode array float dd, dr; // translation from first to second diode array. int a,b,c; // renormalize dra, ddec to have a norm of 1 arc min on sky { double nrm; nrm = dra*15*cos(dec * 3.1415926/180); nrm = sqrt(ddec*ddec + nrm*nrm)*60; ddec /= nrm; dra /= nrm; } // orthogonal, vertical, vector. up-down depends on spin orientation. //ora = ddec / (15*cos(dec * 3.1415926/180)); //odec = - dra * (15*cos(dec * 3.1415926/180)); // size of diode array // along spin, 4*1.88 arc minutes a1 = 4*1.88*dra; d1 = 4*1.88*ddec; // orthogonal, 46 diodes of 1.88 arc minutes each // (mail from Silva Masi, 1.41 deg in elevation between center of pixel 1 // and center of pixel 46). a2 = (1.88*46)*ora; d2 = (1.88*46)*odec; alow = ra - a1/2 - a2/2; dlow = dec - d1/2 - d2/2; // shift with second diode array dd = 10. * ddec - .5 * odec; dr = 10. * dra - .5 * ora; // Fin min/max values of delta, and thus imin imax in gsc.. dmin = dmax = dlow; if ((dlow+d1) < dmin) dmin = dlow+d1; if ((dlow+d2) < dmin) dmin = dlow+d2; if ((dlow+d1+d2) < dmin) dmin = dlow+d1+d2; if ((dlow+d1) > dmax) dmax = dlow+d1; if ((dlow+d2) > dmax) dmax = dlow+d2; if ((dlow+d1+d2) > dmax) dmax = dlow+d1+d2; if ((dmin+dd) < dmin) dmin = dmin+dd; if ((dmax+dd) > dmax) dmax = dmax+dd; a=0; c=nstars-1; while (a+1=0) diodes[(int)(f*46)] = 1000-stars[c].mag; if (f>=0) { int idiod = (int)(f*46); double flux = 3500*pow(10,-(stars[c].mag/250.)); if (secondes - laststartime[idiod][9]>dt/3 || laststarmag[idiod][9] != flux) { // mag=0 -> flux = 3500 (Vega) // decalage... for (a=0; a<9; a++) { laststartime[idiod][a] = laststartime[idiod][a+1]; laststarmag[idiod][a] = laststarmag[idiod][a+1]; } laststartime[idiod][9] = secondes; laststarmag[idiod][9] = flux; } } // Test dans la deuxieme barrette if (sst2Bars) { f = isStarInArray(stars[c].ra, stars[c].dec, alow+dr, dlow+dd, a1, d1, a2, d2); //if (f>=0) diodes[(int)(f*46)] = 1000-stars[c].mag; if (f>=0) { int idiod = (int)(f*46); double flux = 3500*pow(10,-(stars[c].mag/250.)); if (secondes - laststartime[idiod][9]>dt/3 || laststarmag[idiod][9] != flux) { for (a=0; a<9; a++) { laststartime[idiod][a] = laststartime[idiod][a+1]; laststarmag[idiod][a] = laststarmag[idiod][a+1]; } laststartime[idiod][9] = secondes; laststarmag[idiod][9] = flux; } } } } for (c=0; c<46; c++) { double signal = 0; for (a=0; a<10; a++) signal += sig2(secondes - laststartime[c][a]) * laststarmag[c][a]; signal += (c*2-46); // offset specifique a chaque diode if (signal > 2047) signal = 2047; if (signal < -2047) signal = -2047; //if (signal < 0) signal += 4096; signal += 2048; diodes[c] = signal; //if ( (secondes - laststartime[c])
0) diodes[c-1] = laststarmag[c]; //diodes[c] = (secondes - laststartime[c])
>4; t = mode_transmission_telemesure[t][type]; //e = mode_transmission_flash[e][type]; a=0; if( (t >0) && ( (numblock%t)==0) ) a=1; //if( (e >0) && ( (numblock%e)==0) ) a+=2; return(a); } int SimulMissionBloc(tmtc* tt) { double t; static int iLast = -99; if (tickStart<0) { restartTime(); SimBlocParam(tt); return 1; } if (iLast < 0) { iLast = 0; SimBlocReglage(tt); return 1; } //if (stateChanged) { // emission d'un nouveau bloc reglage // SimBlocReglage(tt); // stateChanged=0; // return 1; //} // Avons-nous un bloc bolo a transmettre ? // duree contenue dans un bloc * numero de bloc... t = numblock * per_echant() * (double)nb_per_block*2.; if (iLast == 0 && ((TickCount() - tickStart)/60. < t/overSpeed)) return 0; if ((TickCount() - tickStart)/60. > t/overSpeed + 10) { numblock = (TickCount() - tickStart)/60.*overSpeed / (per_echant() * (double)nb_per_block*2.); iLast = 0; } if (iLast == 0) { numblock++; // On construit un scan... PerformScan(); } // On decide des blocs a envoyer... //#define need(n) \ // (((simReglage.vitesse[n])>0) && ((numblock % (simReglage.vitesse[n])) ==0)) //#define need(n) \ // (((vitesse[n])>0) && ((numblock % (vitesse[n])) ==0)) #define need(n) \ ((emission_blocks(n))!=0) if (iLast < 1 && need(block_param)) { iLast = 1; SimBlocParam(tt); return 1; } if (iLast < 2 && need(block_reglage)) { iLast = 2; SimBlocReglage(tt); return 1; } if (iLast < 4 && need(block_gps)) { iLast = 4; SimBlocGPS(tt); return 1; } if (iLast < 5 && need(block_une_periode)) { iLast = 5; SimBloc1Per(tt); return 1; } if (iLast < 6 && need(block_bolo)) { iLast = 6; SimBlocBolo(tt); return 1; } if (iLast < 7 && need(block_bolo_comprime)) { iLast = 7; SimBlocBoloComp(tt); return 1; } if (iLast < 8 && need(block_sst)) { iLast = 8; SimBlocSST(tt); return 1; } if (iLast < 9 && need(block_sst_comprime)) { iLast = 9; SimBlocSSTComp(tt); return 1; } iLast = 0; return 0; } #define verif_bolo {bolo=tc[2];if(bolo>nb_max_bolo) {printf("erreur: telecommande pour bolo=%d >nb_max_bolo \n",bolo);return;}} void simul_telecommandes_reduites(unsigned char* tc) { int bolo,a;//,i; if(tc[0]!=tc_reduite) return; stateChanged = 1; switch(tc[1]) { case tc2_bolo_dacV: verif_bolo;a=dac_V(simReglage.bolo[bolo]); a=new_val_dac(a,tc[3]); simReglage.bolo[bolo].mot1=(simReglage.bolo[bolo].mot1&0xfff000ff) | (a<<8); printf("Simul: nouvelle valeur dacV=%d \n",a); break; case tc2_bolo_dacI: verif_bolo;a=dac_I(simReglage.bolo[bolo]); a=new_val_dac(a,tc[3]); simReglage.bolo[bolo].mot1=(simReglage.bolo[bolo].mot1&0x000fffff) | (a<<20); printf("Simul: nouvelle valeur dacI=%d \n",a); break; case tc2_bolo_gain: verif_bolo; simReglage.bolo[bolo].mot1=(simReglage.bolo[bolo].mot1&0xffffff00) | tc[3]; printf("Simul: nouvelle valeur gain=%x \n",tc[3]); break; case tc2_bolo_dacT: verif_bolo;a=dac_T(simReglage.bolo[bolo]); a=new_val_dac(a,tc[3]); simReglage.bolo[bolo].mot2=(simReglage.bolo[bolo].mot2&0xfff000ff) | (a<<8); printf("Simul: nouvelle valeur dacT=%d \n",a); break; case tc2_bolo_dacL: verif_bolo;a=dac_L(simReglage.bolo[bolo]); a=new_val_dac(a,tc[3]); simReglage.bolo[bolo].mot2=(simReglage.bolo[bolo].mot2&0x000fffff) | (a<<20); printf("Simul: nouvelle valeur dacL=%d \n",a); break; case tc2_bolo_voie: verif_bolo; simReglage.bolo[bolo].mot2=(simReglage.bolo[bolo].mot2&0xffffff00) | tc[3]; printf("Simul: nouvelle valeur voie=%d \n",tc[3]); break; case tc2_horloge: switch(tc[2]) { case tc3_commande : printf(" commande generale numero =%d \n",tc[3]); //commande_directe(tc[3]); break; case tc3_periode : simReglage.horloge.periode=tc[3]; printf("Simul: nouvelle valeur periode=%d \n",tc[3]); break; case tc3_nb_mesures : simReglage.horloge.nb_mesures=tc[3]; printf("Simul: nouvelle valeur nb_mesures=%d \n",tc[3]); break; case tc3_temp_mort : simReglage.horloge.temp_mort=tc[3]; printf("Simul: nouvelle valeur temp_mort=%d \n",tc[3]); break; case tc3_flag : simReglage.horloge.flag=tc[3]; printf("Simul: nouvelle valeur flag=%d \n",tc[3]); break; default : // i=tc[2]-tc3_vitesse; // if( (i>=0) && (i codage def_gains #define gainbrut(aa) ((char)(((aa).mot1&0x1f))) #define gain(aa) gains_reels[gainbrut(aa)] #define codsignal(muv) (muv/1.e5 * (65536.*gain(simReglage.bolo[i]))) void PerformScan() { int nb_coups,aa; int i,j; float ra,dec,dra,ddec; struct localgeom lastGeom; struct localgeom nextGeom; double perech; double secondes; nb_coups= simReglage.horloge.nb_mesures/2 - simReglage.horloge.temp_mort; aa = (nb_coups<<14) + (nb_coups*190) ; perech = per_echant(); for (j=0; j> 1 )+aa; } calcLocalGeom(&lastData, &lastGeom); calcLocalGeom(&nextData, &nextGeom); dra = nextGeom.ra-lastGeom.ra; ddec = nextGeom.dec-lastGeom.dec; // renormalize dra, ddec to have a norm of 1 arc min on sky lastGeom.dra = dra; lastGeom.ddec = ddec; if (lastGeom.dra>12) lastGeom.dra-=24; if (lastGeom.dra<-12) lastGeom.dra+=24; { double nrm; nrm = lastGeom.dra*15*cos(lastGeom.dec * 3.1415926/180); nrm = sqrt(lastGeom.ddec*lastGeom.ddec + nrm*nrm)*60; lastGeom.ddec /= nrm; lastGeom.dra /= nrm; } ra = lastGeom.ra + dra*iGensample/oversample; dec = lastGeom.dec + ddec*iGensample/oversample; // if (iGensample == 25 && iBuf == 0 && !hasDebugInfo) { // hasDebugInfo=ra*1000; // } /* if (fabs(ra-17.41980)<0.00002 && fabs(dec-85.95694)<0.00002) { //strcpy(debugInfo, " Debug 1"); hasDebugInfo=iGensample; }*/ SSTSignal(ra,dec,lastGeom.dra,lastGeom.ddec,lastGeom.ora,lastGeom.odec,datadiod[j],secondes); latitude = lastData.pos.lat + (nextData.pos.lat-lastData.pos.lat)*iGensample/oversample; longitude = lastData.pos.lon + (nextData.pos.lon-lastData.pos.lon)*iGensample/oversample; mjd = lastData.mjd + (nextData.mjd-lastData.mjd)*iGensample/oversample; //JD+YGH on change le tps entre 2 echantillons pour simuler une periode differente //mjd = lastData.mjd + (nextData.mjd-lastData.mjd)*iGensample/oversample/1.5; iGensample ++; if (iGensample >= oversample) { iGensample = 0; iBuf++; if (iBuf >= szmisdatbuf-1) { needBuf = curBuf; curBuf = 1-curBuf; iBuf = 0; } } } } void SimBlocParam(tmtc* tt) { /* simParam.n_max_bolo = nb_max_bolo; simParam.n_per_block = nb_per_block; simParam.n_max_mes_per = nb_max_mes_per; simParam.nb_bolo = 6; for (i=0; ivi.btt))->param = simParam; valide_block(&tt->vi.btt,block_param,numblock); } void SimBlocReglage(tmtc* tt) { ((block_type_reglage*)(&tt->vi.btt))->reglage = simReglage; valide_block(&tt->vi.btt,block_reglage,numblock); } void SimBlocBolo(tmtc* tt) { int i,j; block_type_bolo* blk = (block_type_bolo*)(&tt->vi.btt); for (j=0; jdata_bolo[i][j] = databolo[i][j]; } } for(i=5;idata_bolo[i][j] = 0; } valide_block(&tt->vi.btt,block_bolo,numblock); } void SimBlocBoloComp(tmtc* tt) { block_type_bolo_comprime blk2; block_type_bolo* blk; int i; SimBlocBolo(tt); blk = (block_type_bolo*) &(tt->vi.btt); for(i=0;idata_bolo[i],(unsigned long*)blk2.data_bolo[i],nb_per_block*2,1); } valide_block((block_type_modele*)&blk2,block_bolo_comprime,numero_block(blk)); memcpy(&tt->vi.btt,&blk2,sizeof(blk2)); } void SimBloc1Per(tmtc* tt) { int i,j; block_type_une_periode* blk = (block_type_une_periode*)(&tt->vi.btt); for (i=0; i<6; i++) for (j=0; jbol_per[i][j] = (value + 0x4000) ^ 0x7fff; } for (i=6; ibol_per[i][j] = 0; valide_block(&tt->vi.btt,block_une_periode,numblock); } // Format bloc GPS // $GPGGA,hhmmss.ss,ddmm.mmmm,n,dddmm.mmmm,e,q,ss,y.y,a.a,z, void SimBlocGPS(tmtc* tt) { int h,m,s; int dlo, dla; double mlo, mla, j; char clo, cla; block_type_gps* blk = (block_type_gps*)(&tt->vi.btt); j = (mjd - (int)mjd + 0.5) * 24; h = j; j = (j-h)*60; m = j; s = (j-m)*60; cla = 'N'; j = latitude; if (latitude < 0) { cla = 'S'; j=-latitude; } dla = j; mla = (j-dla)*60; clo = 'E'; j = longitude; if (longitude < 0) { clo = 'W'; j=-longitude; } dlo = j; mlo = (j-dlo)*60; //strcpy(blk->gps, "$GPGGA,042232,3827.7653,N,00134.2222,E,1,07,2.3,32310.3,M\n"); //sprintf(blk->gps, "$GPGGA,%02d%02d%02d,%02d%07.4f,%1c,%02d%07.4f,%1c,1,07,2.3,32310.3,M\n", sprintf(blk->gps, "$%02d%02d%02d,%02d%07.4f,%1c,%02d%07.4f,%1c,1,07,02.3,000,32310.3,M,32280,M,,\n", h,m,s, dla, mla, cla, dlo, mlo, clo); valide_block(&tt->vi.btt,block_gps,numblock); } void code_sst(block_type_sst* blk, int i, int* diodes); // diodes = tableau a 48 entrees // diodpermut[i] = channel de la diode i static int diodpermut[46]= { 8,24,40, 9,25,41,10,26,42,11, 27,43,16,32, 1,17,33, 2,18,34, 3,19,35,12,28,44,13,29,45,14, 30,46,15,31,47,20,36, 5,21,37, 6,22,38, 7,23,39}; // voies 0 et 4 non connectees, voie 1 en panne. void code_sst(block_type_sst* blk, int i, int* diodes) { int j; // 0-5 : numero du bloc de 8 diodes int k; // 0-2 : indice du bloc de 4 bits (une diode = 12 bits = 3 blocs de 4 bits), MSB=0 int l; // 0-7 : indice de la diode dans son bloc (8 diodes * 4 bits = 1 mot de 32 bits) int dd[46]; // numero de la diode (0-47) = j*8+l; // indice dans le bloc sst du mot de 32 bits (0-17) = j*3+k; // indice dans mot de 32 bits du premier bit utile = 4*l; // Permutation des diodes for (j=0; j<46; j++) { dd[j] = diodes[j]; diodes[j]=0; } diodes[46] = diodes[47] = 0; for (j=0; j<46; j++) { diodes[diodpermut[j]] = dd[j]; } diodes[0] = 23 + 2048; // canal non connecte diodes[4] = 45 + 2048; // canal non connecte diodes[1] = 1234 + 2048; // canal 1 en panne pour le moment. for (j=0; j<6; j++) for (k=0; k<3; k++) { blk->sst[i][j*3+k] = 0; for (l=0; l<8; l++) { short bit4 = (diodes[j*8+l] >> 4*(2-k)) & 0xF; blk->sst[i][j*3+k] += (bit4 << (4*l)); } } } void SimBlocSST(tmtc* tt) { int i; block_type_sst* blk = (block_type_sst*)(&tt->vi.btt); for (i=0; ivi.btt,block_sst,numblock); } #define place_paquet(i,j) ((i/8) * 24 + j*8 + (i%8) ) void SimBlocSSTComp(tmtc* tt) { block_type_sst_comprime blk2; block_type_sst* blk; int j,k,jc; unsigned long sst_vrai[nb_per_block*2]; unsigned long a,b0,b1,b2; SimBlocSST(tt); blk = (block_type_sst*) &(tt->vi.btt); /* pour nb_per_block=36 periodes completes on a 72 points a comprimer pour chaque diode */ jc=0; for(j=0;j<48;j++) /* jc = bolo_comprime j=bolo normal */ { if( (j!=0) && (j!=4) ) { for(k=0;ksst[k][a/8] >>( (a%8)*4) ) & 0xf; a=place_paquet(j,1); b1= ( blk->sst[k][a/8] >>( (a%8)*4) ) & 0xf; a=place_paquet(j,2); b2= ( blk->sst[k][a/8] >>( (a%8)*4) ) & 0xf; sst_vrai[k]=( (b0<<8) | (b1<<4) | b2 ) ; } compress_4_1(sst_vrai,blk2.sst[jc],nb_per_block*2,1); jc++; } } valide_block((block_type_modele*)&blk2,block_sst_comprime,numero_block(blk)); memcpy(&tt->vi.btt,&blk2,sizeof(blk2)); } #define SIDRATE 0.9972695677 #define refTS (13*3600. + 10.*60. + 33.1) #include "aa_hadec.h" void calcLocalGeom(struct misdat* md, struct localgeom* geom) { // Calcul temps sideral local // GMT TS Pour JD = 2450000 = 13h10m33.1s double ts; // secondes double ha, dec, lat, alt, az; // La longitude est vers l'EST ts = md->mjd * 86400. / SIDRATE + refTS + (md->pos.lon/15. * 3600.); ts = ts - 86400.*(int)(ts/86400); // Elevation = 41 degres... -> Ra & Dec alt = 41. * M_PI/180; lat = md->pos.lat * M_PI/180; az = md->pos.azimut * M_PI/180; aa_hadec(lat, alt, az, &ha, &dec); geom->ra = - (ha * 180. / M_PI / 15) + (ts/3600.); geom->dec = (dec * 180. / M_PI); alt = (41. + 1./60.) * M_PI/180; aa_hadec(lat, alt, az, &ha, &dec); geom->ora = - (ha * 180. / M_PI / 15) + (ts/3600.) - geom->ra; geom->odec = (dec * 180. / M_PI) - geom->dec; // en fait, prendre plus tard l'orthogonal du vecteur de balayage, dans le // sens le plus proche de cette valeur (indetermination sens de rotation) // pour tenir compte de la pendulation } // Etalonnage approx : 3e8 V/W // Cte temps 10 ms #define DeuxPi (2* M_PI) // Generation aleatoire gaussienne de centre "am" et de sigma "s". double GauRnd(double am, double s) { double x,A,B; LAB10: A = drand01(); if ( A == 0. ) goto LAB10; B = drand01(); x = am + s * sqrt(-2.*log(A))*cos(DeuxPi*B); return(x); }