source: Sophya/trunk/AddOn/TAcq/svv2mtx.cc@ 3699

Last change on this file since 3699 was 3699, checked in by ansari, 16 years ago

recherche du probleme nb de paquets lus=50% nb de paquets, Reza 4/12/09

File size: 7.1 KB
Line 
1// Utilisation de SOPHYA pour faciliter les tests ...
2#include "sopnamsp.h"
3#include "machdefs.h"
4
5/* ----------------------------------------------------------
6 Programme de lecture multi canaux pour BAORadio
7 R. Ansari, C. Magneville
8 V : Mai 2009
9 ---------------------------------------------------------- */
10
11// include standard c/c++
12#include <math.h>
13#include <stdio.h>
14#include <stdlib.h>
15#include <string.h>
16
17#include <iostream>
18#include <string>
19
20#include "pexceptions.h"
21#include "tvector.h"
22#include "fioarr.h"
23// #include "tarrinit.h"
24#include "ntuple.h"
25#include "histinit.h"
26#include "matharr.h"
27#include "timestamp.h"
28
29// include sophya mesure ressource CPU/memoire ...
30#include "resusage.h"
31#include "ctimer.h"
32#include "timing.h"
33
34
35int Usage(bool fgshort=true);
36// Pour traitement (calcul FFT et visibilites (ProcA) 1 fibre, 2 voies RAW)
37int ProcSVFilesVJun09(string& inoutpath, int imin, int imax, int istep, int jf1, int jf2, int nfreq, int card=1);
38int ProcSVFiles(string& inoutpath, int imin, int imax, int istep, int jf1, int jf2, int nfreq);
39
40//----------------------------------------------------
41//----------------------------------------------------
42int main(int narg, char* arg[])
43{
44 if ((narg>1)&&(strcmp(arg[1],"-h")==0)) return Usage(false);
45 if (narg<4) return Usage(true);
46
47 HiStatsInitiator _inia;
48 // TArrayInitiator _inia;
49
50 int rc = 0;
51 try {
52 string inoutpath = arg[1];
53 int imin=0;
54 int imax=0;
55 int istep=1;
56 sscanf(arg[2],"%d,%d,%d",&imin,&imax,&istep);
57 int jf1=0;
58 int jf2=0;
59 int nfreq=0;
60 sscanf(arg[3],"%d,%d,%d",&jf1,&jf2,&nfreq);
61 int card=1;
62 if (narg>4) card=atoi(arg[4]);
63 cout << " ----- svv2mtx.cc Start - InOutPath= " << inoutpath << " IMin,Max,Step="
64 << imin << "," << imax << "," << istep << " Card=" << card << endl;
65 cout << "Frequency num range JF=" << jf1 << "," << jf2 << "," << nfreq << " ------- " << endl;
66 ResourceUsage resu;
67 rc=ProcSVFilesVJun09(inoutpath, imin, imax, istep, jf1, jf2, nfreq, card);
68 resu.Update();
69 cout << resu;
70 }
71 catch (PException& exc) {
72 cerr << " svv2mtx.cc catched MiniFITSException " << exc.Msg() << endl;
73 rc = 77;
74 }
75 catch (std::exception& sex) {
76 cerr << "\n svv2mtx.cc std::exception :"
77 << (string)typeid(sex).name() << "\n msg= "
78 << sex.what() << endl;
79 rc = 78;
80 }
81 catch (...) {
82 cerr << " svv2mtx.cc catched unknown (...) exception " << endl;
83 rc = 79;
84 }
85
86 cout << ">>>> svv2mtx.cc ------- END ----------- RC=" << rc << endl;
87 return rc;
88
89}
90
91
92
93// Pour traitement (calcul FFT et visibilites (ProcA) 1 fibre, 2 voies RAW)
94int ProcSVFilesVJun09(string& inoutpath, int imin, int imax, int istep, int jf1, int jf2, int nfreq, int card)
95{
96 Timer tm("ProcSVFiles");
97 char fname[512];
98// NTuple
99 const char* nnames[10] = {"fcsm","ttsm","jfreq","s1","s2","s12","s12re","s12im","s12phi","s12mod"};
100 NTuple nt(10, nnames);
101 double xnt[15];
102 uint_4 nmnt = 0;
103 double ms1,ms2,ms12,ms12re,ms12im,ms12phi,ms12mod;
104
105 TMatrix<r_4> s1, s2;
106 TMatrix<r_4> v12re, v12im, v12phi,v12mod;
107 sa_size_t ncols = (imax-imin+1)/istep;
108 sa_size_t nrows = 10;
109 sa_size_t kc=0;
110 for(int ifile=imin; ifile<=imax; ifile+=istep) {
111 if (card==2)
112 sprintf(fname, "%s/Ch34_%d.ppf",inoutpath.c_str(),ifile);
113 else
114 sprintf(fname, "%s/Ch12_%d.ppf",inoutpath.c_str(),ifile);
115 cout << " ProcSVFiles[" << ifile << "] opening file " << fname << endl;
116 PInPersist pin(fname);
117 string tag1="specV1";
118 string tag2="specV2";
119 string tag12="visiV12";
120 if (card==2) {
121 tag1 = "specV3";
122 tag2 = "specV4";
123 tag12="visiV34";
124 }
125 TVector<r_4> sv1;
126 TVector<r_4> sv2;
127 TVector< complex<r_4> > vv12;
128 pin >> PPFNameTag(tag1) >> sv1;
129 pin >> PPFNameTag(tag2) >> sv2;
130 pin >> PPFNameTag(tag12) >> vv12;
131 if (ifile==imin) {
132 nrows = sv1.Size();
133 cout << " ProcSVFilesVJun09/Info: Output s1,s2 matrix size NRows=NFreq="
134 << nrows << " NCols=NFiles=" << ncols << endl;
135 s1.SetSize(nrows, ncols);
136 s2.SetSize(nrows, ncols);
137 nrows = vv12.Size();
138 cout << " ProcSVFilesVJun09/Info: Output v12 matrix size NRows=NFreq="
139 << nrows << " NCols=NFiles=" << ncols << endl;
140 v12re.SetSize(nrows, ncols);
141 v12im.SetSize(nrows, ncols);
142 v12phi.SetSize(nrows, ncols);
143 v12mod.SetSize(nrows, ncols);
144 }
145 s1.Column(kc) = sv1;
146 s2.Column(kc) = sv2;
147 v12re.Column(kc) = real(vv12);
148 v12im.Column(kc) = imag(vv12);
149 v12phi.Column(kc) = phase(vv12);
150 v12mod.Column(kc) = module(vv12);
151
152// Calcul moyenne dans des bandes en frequence
153 int deltajf=(jf2-jf1)/nfreq;
154 if (deltajf<1) deltajf=1;
155 for(int kf=0; kf<nfreq; kf++) {
156 sa_size_t jfstart=jf1+kf*deltajf;
157 sa_size_t jfend=jfstart+deltajf;
158 if (jfend>jf2) break;
159 nmnt=0; ms1=ms2=ms12=ms12re=ms12im=ms12phi=ms12mod=0.;
160 for(sa_size_t jf=jfstart; jf<jfend; jf++) {
161 ms1 += s1(jf,kc);
162 ms2 += s2(jf,kc);
163 ms12re += v12re(jf,kc);
164 ms12im += v12im(jf,kc);
165 ms12phi += v12phi(jf,kc);
166 ms12mod += v12mod(jf,kc);
167 }
168 nmnt = (jfend-jfstart);
169 if (nmnt>0) {
170 double fnorm = (double)nmnt;
171 xnt[0] = ((int_8)(sv1.Info()["StartFC"])+(int_8)(sv1.Info()["EndFC"]))*0.5;
172 xnt[1] = ((int_8)(sv1.Info()["StartTT"])+(int_8)(sv1.Info()["EndTT"]))*0.5;
173 xnt[2] = kf;
174 xnt[3] = ms1/fnorm;
175 xnt[4] = ms2/fnorm;
176 xnt[5] = ms12/fnorm;
177 xnt[6] = ms12re/fnorm;
178 xnt[7] = ms12im/fnorm;
179 xnt[8] = ms12phi/fnorm;
180 xnt[9] = ms12mod/fnorm;
181 nt.Fill(xnt);
182 }
183 }
184 kc++;
185
186 }
187 if (card==2)
188 sprintf(fname, "%s/Ch34mtx.ppf",inoutpath.c_str());
189 else
190 sprintf(fname, "%s/Ch12mtx.ppf",inoutpath.c_str());
191
192 cout << nt;
193 cout << "ProcSVFilesVJun09: Opening file " << fname << " for writing" << endl;
194 POutPersist po(fname);
195 string tag1="s1";
196 string tag2="s2";
197 string tag12r="v12re";
198 string tag12i="v12im";
199 string tag12p="v12phi";
200 string tagnt="nt12";
201 if (card==2) {
202 tag1="s3";
203 tag2="s4";
204 tag12r="v34re";
205 tag12i="v34im";
206 tag12p="v34phi";
207 tagnt="nt34";
208 }
209 po << PPFNameTag(tag1) << s1;
210 po << PPFNameTag(tag2) << s2;
211 po << PPFNameTag(tag12r) << v12re;
212 po << PPFNameTag(tag12i) << v12im;
213 po << PPFNameTag(tag12p) << v12phi;
214 po << PPFNameTag(tagnt) << nt;
215 cout << "ProcSVFilesVJun09: Matrices s1, s2, v12re, v12im, v12phi, NTuple nt written to file " << fname << endl;
216 return 0;
217}
218
219
220
221
222/* --Fonction-- */
223int Usage(bool fgshort)
224{
225 cout << " --- svv2mtx.cc : Read PPF files produced by mcrd to make matrices BAORadio" << endl;
226 cout << " Usage: mcrd InOutPath Imin,Imax,step NumFreq1,NumFreq2,NBinFreq [card=1]" << endl;
227 if (fgshort) {
228 cout << " mcrd -h for detailed instructions" << endl;
229 return 1;
230 }
231 cout << " InOutPath : Input/Output directory name " << endl;
232 cout << " Imin,Imax,IStep: Input PPF files sequence number \n"
233 << " FileNames=InOutPath/Ch12_II.fits Imin<=II<=Imax II+=IStep \n"
234 << " NumFreq1,NumFreq2,NBinFreq: Freq Zone and number of frequency bins for ntuple\n"
235 << " card=1 Ch12 , card=2 Ch34 " << endl;
236 return 1;
237}
Note: See TracBrowser for help on using the repository browser.