Changeset 2286 for trunk/MagicSoft

Timestamp:

07/17/03 19:02:46 (21 years ago)

Author:

blanch

Message:

Several new features introduced as well as fixed bugs

	- 1/100 events printed out
	- Low gain implemented
	- Different response for outer and inner pixels
	- Some warnings removed
	- pedestal and qe file from inpuit card
	- Faster electronic simulation

File:

• trunk/MagicSoft/Simulation/Detector/Camera/camera.cxx (modified) (45 diffs)

trunk/MagicSoft/Simulation/Detector/Camera/camera.cxx

@@ -21 +21 @@
 //
 // $RCSfile: camera.cxx,v $
-// $Revision: 1.56 $
+// $Revision: 1.57 $
 // $Author: blanch $ 
-// $Date: 2003-06-18 17:13:35 $
+// $Date: 2003-07-17 18:02:46 $
 //
 ////////////////////////////////////////////////////////////////////////
@@ -240 +240 @@
 static float FADC_response_ampl = MFADC_RESPONSE_AMPLITUDE;
 static float FADC_response_fwhm = MFADC_RESPONSE_FWHM;
+static float FADC_resp_ampl_out = MFADC_RESPONSE_AMPLITUDE;
+static float FADC_resp_fwhm_out = MFADC_RESPONSE_FWHM;
 static float FADC_noise = 2.0;
 
@@ -398 +400 @@
   char inname[256];           //@< input file name
   char starfieldname[256];    //@< starfield input file name
+  char qe_filename[256];       //@< qe input file name
 
   char datname[256];          //@< data (ASCII) output file name
@@ -407 +410 @@
 
   char nsbpathname[256];      //@< directory with the dataabse for th NSB
+  char nsbpath_outer[256];    //@< directory with the dataabse for outer pixels
 
   char flag[SIZE_OF_FLAGS + 1];  //@< flags in the .rfl file
-  char flag_new[4];              //@< New flag for the run header in the .rfl file
+  char flag_new[SIZE_OF_FLAGS + 1];              //@< New flag for the run header in the .rfl file
 
   int reflector_file_version;  //@< vrsion of he reflector file
@@ -480 +484 @@
   int flagstoring = 0;
 
-  float fCorrection;          //@< Factor to apply to pixel values (def. 1.)
-
   int ntrigger = 0;           //@< number of triggers in the whole file
   int btrigger = 0;           //@< trigger flag
@@ -493 +495 @@
   float fpixelthres[CAMERA_PIXELS];         //@< Threshold values
 
-  TArrayC *fadcValues;  //@< the analog Fadc signal for pixels
-
-  float degTriggerZone;       //@< trigger area in the camera (radius, in deg.)
+  TArrayC *fadcValues;     //@< the analog Fadc High gain signal for pixels
+  TArrayC *fadcValuesLow;  //@< the analog Fadc Low gain signal for pixels
 
   int still_in_loop = FALSE;
@@ -524 +525 @@
 
   ct_NPixels=camgeom.GetNumPixels();
-  read_QE();
-  read_WC();
 
   //!@' @#### Definition of variables for |getopt()|.
@@ -598 +597 @@
   Data_From_STDIN = get_data_from_stdin();
 
+  // read WC and QE files
+
+  strcpy( qe_filename, get_qe_filename());
+
+  read_QE(qe_filename);
+  read_WC();
+
   // write all images, even those without trigger?
 
@@ -618 +624 @@
   get_secure_threhold(&riseDiskThres,&secureDiskThres);
 
-  get_FADC_properties( &FADC_response_ampl, &FADC_response_fwhm);
+  get_FADC_properties( &FADC_response_ampl, &FADC_response_fwhm,
+                       &FADC_resp_ampl_out, &FADC_resp_fwhm_out);
 
   get_Trigger_properties( &Trigger_gate_length, &Trigger_overlaping_time, &Trigger_response_ampl, &Trigger_response_fwhm);
@@ -639 +646 @@
 
   strcpy( inname, get_input_filename() );
+
   strcpy( starfieldname, get_starfield_filename() );
   strcpy( datname, get_data_filename() );
@@ -645 +653 @@
   strcpy( ct_filename, get_ct_filename() );
   strcpy( nsbpathname, get_nsb_directory() );
+  strcpy( nsbpath_outer, get_nsb_directory_outer() );
 
   // get different parameters of the simulation
@@ -657 +666 @@
   //Parameters for starfield rotation
 
-    get_teles_axis(&zenith, &azimutal); 
     Starfield_rotate = get_starfield_rotate();
     log(SIGNATURE,
@@ -668 +676 @@
 
   log(SIGNATURE,
-      "%s:\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n",
+      "%s:\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n\t%20s:\t%s\n",
       "Filenames",
       "In", inname, 
       "Stars", starfieldname,
-      "NSB database", nsbpathname,
+      "NSB database (inner pixels)", nsbpathname,
+      "NSB database (outer pixels)", nsbpath_outer,
       "CT", ct_filename, 
       "Data", datname,
@@ -790 +799 @@
   Int_t Lev0, Lev1; 
   
-  fadcValues =  new TArrayC(FADC_SLICES);
+  fadcValues    =  new TArrayC(FADC_SLICES);
+  fadcValuesLow =  new TArrayC(FADC_SLICES);
 
   // number of pixels for parameters
@@ -813 +823 @@
                     Trigger_response_fwhm);  //@< A instance of the Class MTrigger 
 
+  // Generate databse for trigger electronic noise
+  Trigger.SetElecNoise(Trigger_noise) ;
+
   // Set Right Discriminator threshold, taking into account trigger pixels
-
-
   Trigger.CheckThreshold(&qThreshold[0]);
 
@@ -857 +868 @@
   }
 
-  MFadc fadc(FADC_response_ampl,FADC_response_fwhm) ; //@< A instance of the Class MFadc
+  MFadc fadc(FADC_response_ampl,FADC_response_fwhm,
+             FADC_resp_ampl_out,FADC_resp_fwhm_out) ; //@< A instance of the Class MFadc
 
   //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
@@ -874 +886 @@
 
   for(i=0;i<CAMERA_PIXELS;i++)
-    input_pedestals[i]=10.5;
-  
+    input_pedestals[i]=get_FADC_pedestal();
+
   fadc.SetPedestals(input_pedestals);
+
+  // Generate databse for the Fadc electronic noise
+  
+  fadc.SetElecNoise(FADC_noise);
 
   // Prepare the raw data output
@@ -987 +1003 @@
 
   //  Fill branches for MMcTrigHeader
-  
+
   if(!Trigger_Loop && Write_McTrig){
 
@@ -1134 +1150 @@
   }
   
-  // set trigger area (deg)
-
-  if ( ! get_trigger_radius( &degTriggerZone ) )
-    degTriggerZone = ( ct_Type == 0 ) ? (5.0) : (5.0);
-
-  if ( ! get_correction( &fCorrection ) )
-    fCorrection = 1.0;
 
   // prepare the NSB simulation
@@ -1146 +1155 @@
   //  Instance of the Mlons class
   MLons lons(Trigger_response_ampl, Trigger_response_fwhm,
-             MFADC_RESPONSE_AMPLITUDE, MFADC_RESPONSE_FWHM);
+             FADC_response_ampl,FADC_response_fwhm);
 
   lons.SetPath(nsbpathname);
+
+  //  Instance of the Mlons class
+  MLons lons_outer(Trigger_response_ampl, Trigger_response_fwhm,
+                   FADC_resp_ampl_out,FADC_resp_fwhm_out);
+
+  lons_outer.SetPath(nsbpath_outer);
 
   if( simulateNSB ){
@@ -1226 +1241 @@
   
   strcpy( flag, "                                        \0" );
+  strcpy( flag_new, "    \0" );
 
   // allocate space for PMTs numbers of pixels
@@ -1264 +1280 @@
 
         //  We break the main loop
-        cout<<"Warning: Expected start of run header flag, but found:"<<flag<<endl;
+        cout<<"Warning: Expected start of run header flag, but found:"<<flag_new<<endl;
         cout<<"         We break the main loop"<<endl;
         break;      
@@ -1277 +1293 @@
 
       fread( flag, SIZE_OF_FLAGS, 1, inputfile );
-
       while( isA( flag, FLAG_START_OF_EVENT   )){ // while there is a next event
         fread( flag_new, 4, 1, inputfile );
@@ -1284 +1299 @@
 
         //  We break while events loop
-          cout<<"Warning: Expected start of event header flag, but found:"<<flag<<endl;
+          cout<<"Warning: Expected start of event header flag, but found:"<<flag_new<<endl;
           cout<<"         We break the while events loop"<<endl;
           break;      
@@ -1512 +1527 @@
           for(UInt_t ui=0;ui<camgeom.GetNumPixels();ui++)
             if(nsb_phepns_rotated[ui]>0.0){
-              k=lons.GetResponse(nsb_phepns_rotated[ui],0.1,
-                                 & nsb_trigresp[0],& nsb_fadcresp[0]);
+              if(camgeom[ui].GetR()>camgeom[0].GetR()){
+                k=lons_outer.GetResponse(nsb_phepns_rotated[ui],0.1,
+                                   & nsb_trigresp[0],& nsb_fadcresp[0]);
+              }
+              else{
+                k=lons.GetResponse(nsb_phepns_rotated[ui],0.1,
+                                   & nsb_trigresp[0],& nsb_fadcresp[0]);
+              }
               if(k==0){
                 cout << "Exiting.\n";
@@ -1528 +1549 @@
         }
 
+        ntcph+=ncph;
         if ((nshow+ntshow)%100 == 1){
           log(SIGNATURE, "End of this event: %d cphs(+%d). . .\n",
@@ -1535 +1557 @@
           cout<<inumphensb<<endl;
         }
-        ntcph += ncph;
 
         // skip it ?
@@ -1567 +1588 @@
         // TRIGGER HERE
 
-        //
-        //   now the noise of the electronic 
-        //   (preamps, optical transmission,..)  is introduced. 
-        //   This is done inside the class MTrigger by the method ElecNoise. 
-        //   
-        if (addElecNoise){      
-            Trigger.ElecNoise(Trigger_noise) ;
-            
-            fadc.ElecNoise(FADC_noise) ;
-        }
-
         //  We should simulate the AC coupling behaviour:
         //  For the FADC it is only done for the NSB while producinr
@@ -1583 +1593 @@
         //  For the trigger is done in the Trigger.Diskriminate(), which
         //  is called later (it should be seperated to speed up the program
+        //
+
+        //   now the noise of the electronic 
+        //   (preamps, optical transmission,..)  is introduced. 
+        //   This is done inside the class MTrigger by the method ElecNoise. 
+        //   
+
+        if (addElecNoise && nphe2NSB<=inumphe){ 
+          Trigger.ElecNoise(Trigger_noise) ;
+            
+          fadc.ElecNoise(FADC_noise) ;
+        }
 
         //  now a shift in the fadc signal due to the pedestlas is 
@@ -1693 +1715 @@
                           for (j=0;j<FADC_SLICES;j++){
                             fadcValues->AddAt(fadc.GetFadcSignal(i,j),j);
+                            fadcValuesLow->AddAt(fadc.GetFadcLowGainSignal(i,j),j);
                           }
                           EvtData[iconcount]->AddPixel(i,fadcValues,0);
+                          EvtData[iconcount]->AddPixel(i,fadcValuesLow,kTRUE);
                         }
                       }
@@ -1725 +1749 @@
               if (reflector_file_version<6){
                 mcevth.get_times(&ftime, &ltime);
-                McEvt->Fill( -1.0, 
+                McEvt->Fill( 0, 
                              (UShort_t) mcevth.get_primary() , 
                              mcevth.get_energy(),
@@ -1750 +1774 @@
                              (UInt_t)mcevth.get_CORSIKA(), 
                              (UInt_t)mcevth.get_AtmAbs(), 
-                             (UInt_t)mcevth.get_MirrAbs()+mcevth.get_OutOfMirr()+mcevth.get_BlackSpot(), 
+                             (UInt_t)(mcevth.get_MirrAbs()+mcevth.get_OutOfMirr()+mcevth.get_BlackSpot()), 
                              (UInt_t) ncph, 
                              (UInt_t) inumphe,
@@ -1762 +1786 @@
                 mcevth_2.get_times(&ftime, &ltime);
                 chi2=mcevth_2.get_NKGfit(&Nmax, &t0, &tmax, &a, &b, &c);
-                McEvt->Fill( mcevth_2.get_evt_number(),
+                McEvt->Fill((UInt_t) mcevth_2.get_evt_number(),
                             (UShort_t) mcevth_2.get_primary() , 
                              mcevth_2.get_energy(),
@@ -1787 +1811 @@
                              (UInt_t)mcevth_2.get_CORSIKA(), 
                              (UInt_t)mcevth_2.get_AtmAbs(), 
-                             (UInt_t)mcevth_2.get_MirrAbs()+mcevth_2.get_OutOfMirr()+mcevth_2.get_BlackSpot(), 
+                             (UInt_t)(mcevth_2.get_MirrAbs()+mcevth_2.get_OutOfMirr()+mcevth_2.get_BlackSpot()), 
                              (UInt_t) ncph, 
                              (UInt_t) inumphe,
@@ -1898 +1922 @@
                   for (j=0;j<FADC_SLICES;j++){
                     fadcValues->AddAt(fadc.GetFadcSignal(i,j),j);
+                    fadcValuesLow->AddAt(fadc.GetFadcLowGainSignal(i,j),j);
                   }
                   EvtData[0]->AddPixel(i,fadcValues,0);
+                  EvtData[0]->AddPixel(i,fadcValuesLow,kTRUE);
                 }
               }
@@ -1911 +1937 @@
               if (reflector_file_version<6){
                 mcevth.get_times(&ftime, &ltime);
-                McEvt->Fill( -1.0, 
+                McEvt->Fill( 0, 
                              (UShort_t) mcevth.get_primary() , 
                              mcevth.get_energy(),
@@ -1936 +1962 @@
                              (UInt_t)mcevth.get_CORSIKA(), 
                              (UInt_t)mcevth.get_AtmAbs(), 
-                             (UInt_t)mcevth.get_MirrAbs()+mcevth.get_OutOfMirr()+mcevth.get_BlackSpot(), 
+                             (UInt_t)(mcevth.get_MirrAbs()+mcevth.get_OutOfMirr()+mcevth.get_BlackSpot()), 
                              (UInt_t) ncph, 
                              (UInt_t) inumphe,
@@ -1948 +1974 @@
                 mcevth_2.get_times(&ftime, &ltime);
                 chi2=mcevth_2.get_NKGfit(&Nmax, &t0, &tmax, &a, &b, &c);
-                McEvt->Fill( mcevth_2.get_evt_number(),
+                McEvt->Fill( (UInt_t) mcevth_2.get_evt_number(),
                              (UShort_t) mcevth_2.get_primary() , 
                              mcevth_2.get_energy(),
@@ -1973 +1999 @@
                              (UInt_t)mcevth_2.get_CORSIKA(), 
                              (UInt_t)mcevth_2.get_AtmAbs(), 
-                             (UInt_t)mcevth_2.get_MirrAbs()+mcevth_2.get_OutOfMirr()+mcevth_2.get_BlackSpot(), 
+                             (UInt_t) (mcevth_2.get_MirrAbs()+mcevth_2.get_OutOfMirr()+mcevth_2.get_BlackSpot()), 
                              (UInt_t) ncph, 
                              (UInt_t) inumphe,
@@ -2088 +2114 @@
         } // end if found end of file
       } // end if found end of run
-
+      
     } // end if else found start of run
   } // end big while loop
@@ -2101 +2127 @@
   get_starfield_center(&sfRaH,&sfRaM,&sfRaS,&sfDeD,&sfDeM,&sfDeS);
   if (reflector_file_version<6){
-    get_teles_axis(&telestheta,&telesphi);
+    telestheta=-10.0;
+    telesphi=-10.0;
     mcevth.get_theta_range(&shthetamin, &shthetamax);
     mcevth.get_phi_range(&shphimin,&shphimax);
@@ -2713 +2740 @@
 //!@{
 void 
-read_QE(void){
+read_QE(char fname[256]){
   ifstream qefile;
   char line[LINE_MAX_LENGTH];
@@ -2724 +2751 @@
   // try to open the file
 
-  log("read_QE", "Opening the file \"%s\" . . .\n", QE_FILE);
-  
-  qefile.open( QE_FILE );
+  log("read_QE", "Opening the file \"%s\" . . .\n", fname);
+  
+  qefile.open( fname );
   
   // if it is wrong or does not exist, exit
   
   if ( qefile.bad() )
-    error( "read_QE", "Cannot open \"%s\". Exiting.\n", QE_FILE );
+    error( "read_QE", "Cannot open \"%s\". Exiting.\n", fname );
   
   // read file
@@ -2846 +2873 @@
   // try to open the file
 
-  log("read_QE", "Opening the file \"%s\" . . .\n", WC_FILE);
+  log("read_WC", "Opening the file \"%s\" . . .\n", WC_FILE);
   
   wcfile.open( WC_FILE );
@@ -3894 +3921 @@
     // store the new photoelectron
 
-    fadc->Fill(ipixnum,(t-*tmin_ns) , trigger->FillShow(ipixnum,t-*tmin_ns));
+    fadc->Fill(ipixnum,(t-*tmin_ns) , 
+               trigger->FillShow(ipixnum,t-*tmin_ns),
+               !((*camgeom)[ipixnum].GetR()>(*camgeom)[0].GetR()));
     
     *itotnphe += 1;