Skip to content

GitLab

Docker-in-Docker (DinD) capabilities of public runners deactivated. More info

Commit 53856325 authored by Reza ANSARI's avatar Reza ANSARI
Browse files

Ajout calcul puissance a 1420 MHz +/- 1 +/- 5 au NTuple de visi2ntac.cc et... 

Ajout calcul puissance a 1420 MHz +/- 1 +/- 5 au NTuple de visi2ntac.cc et trace ds le scipt plac8.pic, Reza 12/03/2015
parent 8c4ec2b2
Hide whitespace changes
Inline Side-by-side
Showing with 113 additions and 9 deletions
plac8.pic 0 → 100644
View file @ 53856325
################## Script plac8.pic #################
## Script pour tracer les 8 auto-correlations de PAON4
if ( $# < 1 ) then
echo 'plac8.pic/missing argument '
echo 'Usage: plac8 input_base_name (without .ppf)'
return
endif
set flnm $1
set titre 'PAON4 - transit'
# Pour les axes
setaxesatt 'font=helvetica,bold,16 fixedfontsize minorticks'
# Ouvrir le fichier
openppf $flnm.ppf
set NT $flnm
print $NT
## Facteur de normalisation entre les voies x 1.e-4
## 1 1/3.48 1.06 0.92 1 0.98 0.377 1.32
set FNORM ( 1.e-4 0.287e-4 1.06e-4 0.92e-4 1.e-4 0.96e-4 0.377e-4 1.32e-4 )
## couleurs pour les traces
set COLS ( black red orange brown blue green skylbue forestgreen )
###### Trace de la puissance large bande
## On trace p1 en fct du temps
newwin 1 1 800 600
plot2d $NT time/3600 p1*$FNORM[0] 1 " nsta cpts notit $COLS[0] "
### trace de p2 ... p8
for i 1:8
j = i+1
plot2d $NT time/3600 p$j*$FNORM[i] 1 "nsta cpts same notit $COLS[i] "
end
### Axes et titres
setaxelabels 'Time (hours)' 'P1...P8 (arb.units)' 'font=helvetica,bolditalic,16 black'
settitle "$titre" ' ' 'font=helvetica,bold,18 black'
newwin 2 2
for i 0:4
j = i+1
plot2d $NT time/3600 p$j*$FNORM[i] 1 "nsta cpts notit $COLS[i] "
settitle "Transit P$j" ' ' 'font=helvetica,bold,14 black'
end
newwin 2 2
for i 4:8
j = i+1
plot2d $NT time/3600 p$j*$FNORM[i] 1 "nsta cpts notit $COLS[i] "
settitle "Transit P$j" ' ' 'font=helvetica,bold,14 black'
end
###### Trace de la puissance a +/-1 MHz @ 1420 MHz
newwin 1 1 800 600
plot2d $NT time/3600 p1HI*$FNORM[0] 1 " nsta cpts notit $COLS[0] "
### trace de p2 ... p8
for i 1:8
j = i+1
plot2d $NT time/3600 p${j}HI*$FNORM[i] 1 "nsta cpts same notit $COLS[i] "
end
### Axes et titres
setaxelabels 'Time (hours)' 'P1...P8 @1420+/-1 MHz ' 'font=helvetica,bolditalic,16 black'
settitle "$titre 1420+/-1 MHz" ' ' 'font=helvetica,bold,18 black'
###### Trace de la puissance a +/-5 MHz @ 1420 MHz
newwin 1 1 800 600
plot2d $NT time/3600 p1HI5*$FNORM[0] 1 " nsta cpts notit $COLS[0] "
### trace de p2 ... p8
for i 1:8
j = i+1
plot2d $NT time/3600 p${j}HI5*$FNORM[i] 1 "nsta cpts same notit $COLS[i] "
end
### Axes et titres
setaxelabels 'Time (hours)' 'P1...P8 @1420+/-5 MHz ' 'font=helvetica,bolditalic,16 black'
settitle "$titre 1420+/-5 MHz" ' ' 'font=helvetica,bold,18 black'
###########################################################
visi2ntac.cc
View file @ 53856325
......@@ -52,7 +52,7 @@ int Usage(bool fgea)
<< " InPathBAO5/6: path for input vismtx_J_iii.ppf on bao5/6 \n"
<< " Imin,Imax : read files for Imin<=iii<=Imax iii+=step \n"
<< " OutPPFFile: Output PPF file name \n"
<< " Jf1,Jf2: frequency index (Jf1<=f<=Jf2) range for power integration, def=768,2768 \n"
<< " Jf1,Jf2: frequency index (Jf1<=f<=Jf2) range for power integration, def=768,2304 \n"
<< " DeltaIAvg: compute average power every DeltaI input vismtx files, def=1 \n"<<endl;
/*
if (fgshort) {
......@@ -76,9 +76,9 @@ int main(int narg, char* arg[])
Istep=1;
string outfile=arg[4];
// frequency range definition
sa_size_t JFmin=768,JFmax=2768;
sa_size_t JFmin=768,JFmax=2304;
int jf1,jf2;
if ((narg>5)&&(strcmp(arg[5],"-")!=0)) sscanf(arg[5],"%d,%d",&jf1,&jf2);
if ((narg>5)&&(strcmp(arg[5],".")!=0)) sscanf(arg[5],"%d,%d",&jf1,&jf2);
JFmin=jf1; JFmax=jf2;
// time averaging interval definition, in term of visibility files
int deltaIavg=1;
......@@ -90,7 +90,13 @@ int main(int narg, char* arg[])
cout << " visi2ntac/Info: Path BAO5:"<<path5<<" BAO6:"<<path6<<"\n"
<< " Imin,max,step="<<Imin<<","<<Imax<<","<<Istep<<" OutFile:"<<outfile<<"\n"
<< JFmin<<" <=NumFreq<= "<<JFmax<<" DeltaIAvg="<<deltaIavg<<" PrtLev="<<prtlev<<endl;
// Bande de frequence de 2 MHz autour du 21 cm (1420 MHz -> JF=2758)
sa_size_t JFmin21=2748,JFmax21=2768; // +/- 1 MHz autour de 1420 MHz
sa_size_t JFmin21_5=2676,JFmax21_5=2840; // +/- 5 MHz autour de 1420 MHz
cout << " pjHI +/-1 MHz @1420 MHz " << JFmin21<<" <=NumFreq<= "<<JFmax21
<< " pjHI +/-5 MHz @1420 MHz " << JFmin21_5<<" <=NumFreq<= "<<JFmax21_5 << endl;
HiStatsInitiator _inia;
// TArrayInitiator _inia;
......@@ -101,12 +107,19 @@ int main(int narg, char* arg[])
paths.push_back(path5);
paths.push_back(path6);
const char* ntnames[10]={"k","time","p1","p2","p3","p4","p5","p6","p7","p8"};
NTuple nt(10,ntnames);
r_8 xnt[20];
const char* ntnames[26]={"k","time",
"p1","p2","p3","p4","p5","p6","p7","p8",
"p1HI","p2HI","p3HI","p4HI","p5HI","p6HI","p7HI","p8HI",
"p1HI5","p2HI5","p3HI5","p4HI5","p5HI5","p6HI5","p7HI5","p8HI5"};
NTuple nt(26,ntnames);
r_8 xnt[30];
// Numero de ligne des auto-correlations dans la matrice des visibilites
sa_size_t KVAC[8]={0,8,15,21,26,30,33,35};
Range freqs(JFmin,JFmax);
Range freqs21(JFmin21,JFmax21);
Range freqs21_5(JFmin21_5,JFmax21_5);
VisiP4Reader vreader(paths, Imin,Imax,Istep,true);
vreader.setPrintLevel(prtlev);
bool fgok=true;
......@@ -122,16 +135,28 @@ int main(int narg, char* arg[])
xnt[0]=cnt; xnt[1]=dateobs.SecondsPart();
if (prtlev>0)
cout << "visi2ntac/Info: dateobs="<<dateobs<<" SecondsPart()="<<dateobs.SecondsPart()<<endl;
for(int k=0; k<8; k++) xnt[2+k]=0.;
for(int k=0; k<24; k++) xnt[2+k]=0.;
}
for(int k=0; k<8; k++) {
for(int k=0; k<8; k++) { // integration en bande large
TVector< complex<r_4> > vac = vismtx.Row(KVAC[k]);
xnt[2+k] += vac(freqs).Sum().real();
}
for(int k=0; k<8; k++) { // integration a +/- 1 MHz @ 1420 MHz
TVector< complex<r_4> > vac = vismtx.Row(KVAC[k]);
xnt[2+8+k] += vac(freqs21).Sum().real();
}
for(int k=0; k<8; k++) { // integration a +/- 5 MHz @ 1420 MHz
TVector< complex<r_4> > vac = vismtx.Row(KVAC[k]);
xnt[2+16+k] += vac(freqs21_5).Sum().real();
}
I++;
if (I==deltaIavg) {
double fnorm=(1./(double)deltaIavg)/(double)(JFmax-JFmin+1);
for(int k=0; k<8; k++) xnt[2+k]*=fnorm;
fnorm=(1./(double)deltaIavg)/(double)(JFmax21-JFmin21+1);
for(int k=0; k<8; k++) xnt[2+8+k]*=fnorm;
fnorm=(1./(double)deltaIavg)/(double)(JFmax21_5-JFmin21_5+1);
for(int k=0; k<8; k++) xnt[2+16+k]*=fnorm;
nt.Fill(xnt);
I=0; cntnt++;
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment