- Timestamp:
- Oct 26, 2011, 11:36:56 AM (13 years ago)
- Location:
- BAORadio/AmasNancay/trunk
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
BAORadio/AmasNancay/trunk/etude_mergeAna.pic
r588 r591 18 18 19 19 echo $nbins1420 20 21 20 22 21 … … 26 25 clearscript IntAt1420 27 26 clearscript sigmaRaw1420Side 27 clearscript sigmaRaw1400a1420 28 28 ############################################################################ 29 29 defscript evolAtCalibFreq … … 150 150 151 151 152 153 154 #col 13 = onoffRaw01420side 155 #col 14 = onoffRaw11420side 152 #col 13 = onoffRaw01420side Ch 0 153 #col 14 = onoffRaw11420side Ch 1 156 154 157 155 ntcol2var onoffevol 13 linCh0All … … 161 159 line2vec vecCh1All $linCh1All 162 160 163 161 #packing values only valid for 500 cycles... 164 162 set packVal "1 10 25 50 100" 165 163 set sigmaCh0 "" … … 188 186 mean1 = ${outvec.sum}/${outvec.size} 189 187 sigma1 = sqrt(${outvec.sumsq}/${outvec.size}-${mean1}*${mean1}) 190 errsig1 = ${sigma 0}/sqrt(2*${outvec.size})188 errsig1 = ${sigma1}/sqrt(2*${outvec.size}) 191 189 set sigmaCh1 "${sigmaCh1} ${sigma1}" 192 190 set errsigCh1 "${errsigCh1} ${errsig1}" … … 222 220 223 221 endscript 222 ################################################## 223 defscript sigmaRaw1400a1420 224 225 c++compile rebining 226 c++link rebining.so dorebin 227 228 229 #col 19 = onoffRaw0f14001420 Ch 0 230 #col 20 = onoffRaw1f14001420 Ch 1 231 232 ntcol2var onoffevol 19 linCh0All 233 line2vec vecCh0All $linCh0All 234 235 ntcol2var onoffevol 20 linCh1All 236 line2vec vecCh1All $linCh1All 237 238 #packing values only valid for 500 cycles... 239 set packVal "1 10 25 50 100" 240 set sigmaCh0 "" 241 set errsigCh0 "" 242 set sigmaCh1 "" 243 set errsigCh1 "" 244 245 246 foreach ipack ( $packVal ) 247 248 #Packing per ipack 249 del invec 250 del outvec 251 cp vecCh0All invec 252 call dorebin $ipack 253 mean0 = ${outvec.sum}/${outvec.size} 254 sigma0 = sqrt(${outvec.sumsq}/${outvec.size}-${mean0}*${mean0}) 255 errsig0 = ${sigma0}/sqrt(2*${outvec.size}) 256 set sigmaCh0 "${sigmaCh0} ${sigma0}" 257 set errsigCh0 "${errsigCh0} ${errsig0}" 258 # 259 del invec 260 del outvec 261 cp vecCh1All invec 262 call dorebin $ipack 263 mean1 = ${outvec.sum}/${outvec.size} 264 sigma1 = sqrt(${outvec.sumsq}/${outvec.size}-${mean1}*${mean1}) 265 errsig1 = ${sigma1}/sqrt(2*${outvec.size}) 266 set sigmaCh1 "${sigmaCh1} ${sigma1}" 267 set errsigCh1 "${errsigCh1} ${errsig1}" 268 end 269 270 271 newnt nt0 x y ey 272 newnt nt1 x y ey 273 274 set x0 ( $packVal ) 275 set y0 ( $sigmaCh0 ) 276 set ey0 ( $errsigCh0 ) 277 278 set x1 ( $packVal ) 279 set y1 ( $sigmaCh1 ) 280 set ey1 ( $errsigCh1 ) 281 282 for i 0:$#x0 283 line2nt nt0 $x0[i] $y0[i] $ey0[i] 284 line2nt nt1 $x1[i] $y1[i] $ey1[i] 285 end 286 287 288 newwin 1 1 289 plot2de nt0 x y 0 ey 1 "blue marker=fcircle,9 notit nsta" 290 plot2de nt1 x y 0 ey 1 "same red marker=fcircle,9 notit nsta" 291 n = ${#x0}-1 292 func $y0[0]/sqrt(x) $x0[0] $x0[n] 100 "same" 293 settitle "Sigma [1400,1420]MHz ${source} Ch 0 (blue) Ch 1 (red)" ' ' $defatt 294 setaxelabels "num of cycles" "Sigma (a.u)" $axedefatt 295 296 297 298 endscript 224 299 ########################################### 300 225 301 defscript evolAt1420SideFreq 226 302 … … 361 437 # Calibration coeff. per cycles 362 438 # 439 440 graphicatt "xylimits=1250,1500,-0.01,0.01" 441 del sonoffovoff0 sonoffovoff1 442 objaoper meanOvOffNoCalib row 0 sonoffovoff0 443 objaoper meanOvOffNoCalib row 1 sonoffovoff1 444 newwin 1 1 445 plot2d sonoffovoff0 (n/8192)*250+1250 val n>0 "blue cpts notit nsta" 446 plot2d sonoffovoff1 (n/8192)*250+1250 val n>0 "red same cpts notit nsta" 447 settitle "Raw (ON-OFF)/OFF ${source} Ch 0 (blue) Ch 1 (red) ${ncycles}cycles" ' ' $defatt 448 setaxelabels "Freq. (MHz)" "I (a.u)" $axedefatt 363 449 364 450 -
BAORadio/AmasNancay/trunk/mergeAnaFiles.cc
r587 r591 111 111 112 112 //------------- 113 //JEC 25/10/11 Perform a median filtering with a sliding window of "halfwidth" half width 114 // It takes care of the edges and is based on the median function (above) 115 void medianFiltering(const TMatrix<r_4> mtx, 116 sa_size_t halfwidth, 117 TMatrix<r_4>& vec) { 118 119 sa_size_t nr = mtx.NRows(); 120 sa_size_t nc = mtx.NCols(); 121 sa_size_t chMin = 0; 122 sa_size_t chMax = nc-1; 123 124 for (sa_size_t ir=0; ir<nr; ir++){ 125 for (sa_size_t ic=0; ic<nc; ic++) { 126 sa_size_t chLow = ic-halfwidth; 127 chLow = (chLow >= chMin) ? chLow : chMin; 128 sa_size_t chHigh = ic+halfwidth; 129 chHigh = ( chHigh <= chMax ) ? chHigh : chMax; 130 TVector<r_4> tmp(mtx(Range(ir),Range(chLow,chHigh)),false); 131 vector<r_4> val; 132 tmp.FillTo(val); 133 vec(ir,ic) = median(val.begin(),val.end()); 134 } 135 } 136 } 137 //------------- 113 138 void split(const string& str, const string& delimiters , vector<string>& tokens) { 114 139 // Skip delimiters at beginning. … … 264 289 265 290 //mean ON-OFF over the list of cycles 291 TMatrix<r_4> meanDiffONOFFovOFF_noCalib(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ); //set to 0 266 292 TMatrix<r_4> meanDiffONOFF_noCalib(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ); //set to 0 267 293 TMatrix<r_4> meanDiffONOFF_perRunCalib(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ); //set to 0 268 294 TMatrix<r_4> meanDiffONOFF_perCycleCalib(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ); //set to 0 269 static const int NINFO= 19;295 static const int NINFO=21; 270 296 char* onffTupleName[NINFO]={"cycle" 271 ,"onoffRaw0","onoffRaw1" 272 ,"onoffRun0","onoffRun1" 273 ,"onoffCycle0","onoffCycle1" 274 ,"onoffRaw01420","onoffRaw11420" 275 ,"onoffRun01420","onoffRun11420" 276 ,"onoffCycle01420","onoffCycle11420" 277 ,"onoffRaw01420side","onoffRaw11420side" 278 ,"onoffRun01420side","onoffRun11420side" 279 ,"onoffCycle01420side","onoffCycle11420side" 297 ,"onoffRaw0","onoffRaw1" 298 ,"onoffRun0","onoffRun1" 299 ,"onoffCycle0","onoffCycle1" 300 ,"onoffRaw01420","onoffRaw11420" 301 ,"onoffRun01420","onoffRun11420" 302 ,"onoffCycle01420","onoffCycle11420" 303 ,"onoffRaw01420side","onoffRaw11420side" 304 ,"onoffRun01420side","onoffRun11420side" 305 ,"onoffCycle01420side","onoffCycle11420side" 306 ,"onoffRaw0f14001420","onoffRaw1f14001420" 280 307 }; 281 308 NTuple onoffevolution(NINFO,onffTupleName); … … 294 321 sa_size_t ch1420bLow = freqToChan(1422); 295 322 sa_size_t ch1420bHigh = freqToChan(1423); 296 323 324 //25/10/11 follow 1400-1420Mhz 325 sa_size_t ch1420 = freqToChan(1420); 326 sa_size_t ch1400 = freqToChan(1400); 297 327 298 328 if (para.debuglev_>0){ … … 631 661 meanDiffONOFF_noCalib += diffOnOff_noCalib; 632 662 663 //JEC 29/10/11 add ON-OFF/OFF 664 TMatrix<r_4> diffOnOffOvOff_noCalib(diffOnOff_noCalib,false); //do not share data 665 TMatrix<r_4> aSpecOffFitltered(NUMBER_OF_CHANNELS,NUMBER_OF_FREQ); 666 sa_size_t halfWidth = 35; //number of freq. bin for the 1/2 width of the filtering window 667 medianFiltering(aSpecOff,halfWidth,aSpecOffFitltered); 668 669 diffOnOffOvOff_noCalib.Div(aSpecOffFitltered); //division in place 670 meanDiffONOFFovOFF_noCalib += diffOnOffOvOff_noCalib; 671 633 672 TMatrix<r_4> diffOnOff_perRunCalib = aSpecOn_BAOCalibRun - aSpecOff_BAOCalibRun; 634 673 meanDiffONOFF_perRunCalib += diffOnOff_perRunCalib; … … 703 742 xnt[18] = (meanInRange_1420aFreq_perCycle(1) + meanInRange_1420bFreq_perCycle(1))/2.; 704 743 744 745 //JEC 25/10/11 follow 1400-1420 MHz bande protege et n'inclue pas le 1420.4Mhz de la Galaxie 746 TVector<r_4> meanInRange_1400a1420Freq_noCalib(NUMBER_OF_CHANNELS); 747 meanInRange(diffOnOff_noCalib,ch1400,ch1420,meanInRange_1400a1420Freq_noCalib); 748 xnt[19] = meanInRange_1400a1420Freq_noCalib(0); 749 xnt[20] = meanInRange_1400a1420Freq_noCalib(1); 750 751 705 752 706 753 //store infos to Ntuple … … 717 764 //Normalisation 718 765 if(totalNumberCycles > 0){ 766 //JEC 29/10 add ON-OFF/OFF 767 meanDiffONOFFovOFF_noCalib /= (r_4)totalNumberCycles; 719 768 meanDiffONOFF_noCalib /= (r_4)totalNumberCycles; 720 769 meanDiffONOFF_perRunCalib /= (r_4)totalNumberCycles; … … 745 794 string tag = "meanNoCalib"; 746 795 fos << PPFNameTag(tag) << meanDiffONOFF_noCalib; 796 797 //JEC 29/10/11 798 tag = "meanOvOffNoCalib"; 799 fos << PPFNameTag(tag) << meanDiffONOFFovOFF_noCalib; 800 747 801 tag = "meanPerRunCalib"; 748 802 fos << PPFNameTag(tag) << meanDiffONOFF_perRunCalib;
Note: See TracChangeset
for help on using the changeset viewer.