// Utilisation de SOPHYA pour faciliter les tests ... #include "sopnamsp.h" #include "machdefs.h" /* ------------------------------------------------------------------ Programme de calcul de spectre moyenne a partir des fichiers fits d'acquisition de BAORadio - donnees brutes ou fft R. Ansari, C. Magneville V1 : Juillet 2008 , V2 : Avril 2009 ------------------------------------------------------------------ */ // include standard c/c++ #include #include #include #include #include #include #include "pexceptions.h" #include "tvector.h" #include "fioarr.h" #include "tarrinit.h" #include "timestamp.h" #include "fftpserver.h" #include "fftwserver.h" #include "FFTW/fftw3.h" // include sophya mesure ressource CPU/memoire ... #include "resusage.h" #include "ctimer.h" #include "timing.h" // include mini-fits lib , et structure paquet BAORadio #include "minifits.h" #include "brpaqu.h" // Definition de type pour les tableaux de float #define DTF r_4 //---- Declaration des fonctions de calcul ---- // Fonction d'analyse 1ere version, pas d'entete ds le fichier, 1 voie int ana_data_0(vector& infiles, string& outfile); // Fonction d'analyse 2eme version , 1 voie / paquet int ana_data_1(vector& infiles, string& oufile, bool fgnotrl=false); // Fonction d'analyse 2eme version , 2 voies / paquet int ana_data_2(vector& infiles, string& oufile, bool fgnotrl=false); // Fonction d'analyse 2eme version , 2 voies / paquet int ana_data_fft_2(vector& infiles, string& oufile, bool fgnotrl=false); //---------------------------------------------------- //---------------------------------------------------- int main(int narg, char* arg[]) { if (narg < 4) { cout << " ---Mean spectrum computation from BAORadio FITS files" << endl; cout << " Usage: mfits2spec ACT OutPPF file1 [file2 file3 ...] " << endl; cout << " Or : mfits2spec ACT OutPPF -infits DirName Imin Imax " << endl; cout << " ACT=-0,-1,-2,-fft2 , -1nt -2nt -fft2nt \n" << " -0: Nancay-July2008 \n" << " -1,-2 : Raw data 1 ou 2 channels / packet(or frame) \n" << " -fft2 : FFT data 2 channels / packet \n" << " nt : (-1nt -2nt -fft2nt) fits files without frame trailer " << endl; cout << " OutPPF : Output PPF file name " << endl; cout << " file1,file2 ... : Input FITS files " << endl; cout << " -infiles DirName Imin Imax : Input fits directory and sequence numbers \n" << " FileNames=DirName/signalII.fits Imin<=II<=Imax " << endl; return 1; } TArrayInitiator _inia; int rc = 0; try { string act = arg[1]; string outppf = arg[2]; vector infiles; if (strcmp(arg[3],"-infits")==0) { if (narg<7) { cout << " mfits2spec/Error arguments - type mfits2spec for help" << endl; return 2; } char nbuff[1024]; char* dirname = arg[4]; int imin = atoi(arg[5]); int imax = atoi(arg[6]); for(int ii=imin; ii<=imax; ii++) { sprintf(nbuff,"%s/signal%d.fits",dirname,ii); infiles.push_back(nbuff); } } else // Liste de noms de fichiers fournis directement for(int i=3; i mai 2009) if ((act=="-1nt")||(act=="-2nt")||(act=="-fft2nt")) fgnotrl=true; // fichier fits SANS trailer if (act == "-0") ana_data_0(infiles, outppf); else if ((act == "-1")||(act=="-1nt")) ana_data_1(infiles, outppf, fgnotrl); else if ((act == "-2")||(act=="-2nt")) ana_data_2(infiles, outppf, fgnotrl); else if ((act == "-fft2")||(act=="-fft2nt")) ana_data_fft_2(infiles, outppf, fgnotrl); else cout << " mfits2spec.cc / Bad argument ACT=" << act << " -> exit" << endl; cout << resu ; } catch (MiniFITSException& exc) { cerr << " mfits2spec.cc catched MiniFITSException " << exc.Msg() << endl; rc = 77; } catch (std::exception& sex) { cerr << "\n mfits2spec.cc std::exception :" << (string)typeid(sex).name() << "\n msg= " << sex.what() << endl; rc = 78; } catch (...) { cerr << " mfits2spec.cc catched unknown (...) exception " << endl; rc = 79; } cout << ">>>> mfits2spec.cc ------- FIN ----------- RC=" << rc << endl; return rc; } inline DTF Zmod2(complex z) { return (z.real()*z.real()+z.imag()*z.imag()); } /*--Nouvelle-Fonction--*/ int ana_data_0(vector& infiles, string& outfile) // Calcul spectre moyen 1 voie, donnees brutes, Donnees Juillet 2008 { TVector spectre; sa_size_t nzm = 0; // Nb de spectres moyennes for(int ifile=0; ifile skipping " << endl; } size_t sx = mff.NAxis1(); size_t sy = mff.NAxis2(); Byte* data = new Byte[sx]; TVector vx(sx); TVector< complex > cfour; FFTPackServer ffts; for(int j=0; j& infiles, string& outfile, bool fgnotrl) // Calcul spectre moyen 1 voie, donnees brutes { TVector spectre; float freq0 = 0; int paqsz = 0; int nfileok=0; sa_size_t nzm = 0; // Nb de spectres moyennes Byte* data = NULL; FFTPackServer ffts; Timer tm("ana_data_1"); size_t totnbytesrd = 0; for(int ifile=0; ifile skipping " << endl; continue; } // Les fichier FITS contiennent l'entet (24 bytes), mais pas le trailer (16 bytes) si fgnotrl=true int incpaqsz=0; if (fgnotrl) { incpaqsz=16; cout << " Warning : FITS files without frame trailers ..." << endl; } if (paqsz == 0) { // premier passage, on fixe la taille de paquet et on alloue le buffer paqsz = mff.NAxis1()+incpaqsz; data = new Byte[paqsz]; for(int ib=0; ib skipping " << endl; continue; } } size_t sx = mff.NAxis1(); size_t sy = mff.NAxis2(); BRPaquet paq(NULL, data, paqsz); TVector vx(paq.DataSize()); TVector< complex > cfour; for(int j=0; j& infiles, string& outfile, bool fgnotrl) // Calcul spectres moyens 2 voies, donnees brutes { TVector specV1, specV2; TVector< complex > cxspecV12; float freq0v1 = 0; float freq0v2 = 0; int paqsz = 0; int nfileok=0; sa_size_t nzm = 0; // Nb de spectres moyennes Byte* data = NULL; FFTPackServer ffts; Timer tm("ana_data_2"); size_t totnbytesrd = 0; for(int ifile=0; ifile skipping " << endl; continue; } // Les fichier FITS contiennent l'entet (24 bytes), mais pas le trailer (16 bytes) si fgnotrl=true int incpaqsz=0; if (fgnotrl) { incpaqsz=16; cout << " Warning : FITS files without frame trailers ..." << endl; } if (paqsz == 0) { // premier passage, on fixe la taille de paquet et on alloue le buffer paqsz = mff.NAxis1()+incpaqsz; cout << " ana_data_2/ Allocating data , PaqSz=" << paqsz << endl; data = new Byte[paqsz]; for(int ib=0; ib skipping " << endl; continue; } } size_t sx = mff.NAxis1(); size_t sy = mff.NAxis2(); BRPaquet paq(NULL, data, paqsz); TVector vx1(paq.DataSize()/2); TVector vx2(paq.DataSize()/2); TVector< complex > cfour1,cfour2; for(int j=0; j((DTF)(nzm),0.); specV1.Info().Comment() = " SpectreMoyenne (Moyenne module^2) - Raw-Voie 1 (/2)"; specV2.Info().Comment() = " SpectreMoyenne (Moyenne module^2) - Raw-Voie 2 (/2)"; specV1.Info()["NMOY"] = specV2.Info()["NMOY"] = nzm; // Nombre de spectres moyennes specV1.Info()["Freq0"] = freq0v1; specV2.Info()["Freq0"] = freq0v2; POutPersist po(outfile); po << PPFNameTag("specV1") << specV1; po << PPFNameTag("specV2") << specV2; po << PPFNameTag("cxspecV12") << cxspecV12; return 0; } inline DTF Zmod2TwoByte(TwoByteComplex z) { return (z.realD()*z.realD()+z.imagD()*z.imagD()); } /*--Nouvelle-Fonction--*/ int ana_data_fft_2(vector& infiles, string& outfile, bool fgnotrl) // Calcul spectres moyens 2 voies, donnees brutes { TVector specV1, specV2; TVector< complex > cxspecV12; float freq0v1 = 0; float freq0v2 = 0; int paqsz = 0; int nfileok=0; sa_size_t nzm = 0; // Nb de spectres moyennes Byte* data = NULL; Timer tm("ana_data_fft_2"); size_t totnbytesrd = 0; for(int ifile=0; ifile skipping " << endl; continue; } // Les fichier FITS contiennent l'entet (24 bytes), mais pas le trailer (16 bytes) si fgnotrl=true int incpaqsz=0; if (fgnotrl) { incpaqsz=16; cout << " Warning : FITS files without frame trailers ..." << endl; } if (paqsz == 0) { // premier passage, on fixe la taille de paquet et on alloue le buffer paqsz = mff.NAxis1()+incpaqsz; cout << " ana_data_fft_2/ Allocating data , PaqSz=" << paqsz << endl; data = new Byte[paqsz]; for(int ib=0; ib skipping " << endl; continue; } } size_t sx = mff.NAxis1(); size_t sy = mff.NAxis2(); BRPaquet paq(NULL, data, paqsz); TVector vx1(paq.DataSize()/2); TVector vx2(paq.DataSize()/2); TVector< complex > cfour1,cfour2; if (!specV1.IsAllocated()) specV1.SetSize(paq.DataSize()/4+1); if (!specV2.IsAllocated()) specV2.SetSize(paq.DataSize()/4+1); if (!cxspecV12.IsAllocated()) cxspecV12.SetSize(paq.DataSize()/4+1); for(int j=0; j((DTF)(nzm),0.); specV1.Info().Comment() = " SpectreMoyenne (Moyenne module^2) - FFT-Voie 1 (/2)"; specV2.Info().Comment() = " SpectreMoyenne (Moyenne module^2) - FFT-Voie 2 (/2)"; specV1.Info()["NMOY"] = specV2.Info()["NMOY"] = nzm; // Nombre de spectres moyennes specV1.Info()["Freq0"] = freq0v1; specV2.Info()["Freq0"] = freq0v2; POutPersist po(outfile); po << PPFNameTag("specV1") << specV1; po << PPFNameTag("specV2") << specV2; po << PPFNameTag("cxspecV12") << cxspecV12; return 0; }