Changeset 1512 for trunk/MagicSoft/Simulation/Detector

Timestamp:

09/04/02 10:57:42 (22 years ago)

Author:

blanch

Message:

Modifications done to use MGeomCam from MARS.

File:

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

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

@@ -21 +21 @@
 //
 // $RCSfile: camera.cxx,v $
-// $Revision: 1.40 $
+// $Revision: 1.41 $
 // $Author: blanch $ 
-// $Date: 2002-07-16 16:15:22 $
+// $Date: 2002-09-04 09:57:42 $
 //
 ////////////////////////////////////////////////////////////////////////
@@ -72 +72 @@
 #include "MMcTrigHeader.hxx"
 #include "MMcFadcHeader.hxx"
+#include "MGeomCamMagic.h"
+#include "MGeomPix.h"
 
 /*!@" 
@@ -495 +497 @@
   UInt_t corsika = 5200 ;
 
-
-  struct camera cam; // structure holding the camera definition
+  MGeomCamMagic camgeom; // structure holding the camera definition
 
 
@@ -759 +760 @@
 
   read_ct_file();
-
-  // read camera setup
-
-  read_pixels(&cam);
 
   Int_t Lev0, Lev1; 
@@ -1112 +1109 @@
     
     k = produce_nsbrates( starfieldname,
-                          &cam,
+                          &camgeom,
                           nsbrate_phepns);
 
@@ -1122 +1119 @@
     // calculate diffuse rate correcting for the pixel size
 
-    for(i=0; i<cam.inumpixels; i++){
-      diffnsb_phepns[i] = meanNSB * 
-        cam.dpixsizefactor[i] * cam.dpixsizefactor[i];
+    const double size_inner = camgeom[0].GetR();
+      
+    for(UInt_t ui=0; ui<camgeom.GetNumPixels(); ui++){
+      const double actual_size = camgeom[ui].GetR();
+      diffnsb_phepns[ui] = meanNSB * 
+        actual_size*actual_size/(size_inner*size_inner);
     }
 
@@ -1134 +1134 @@
         zenfactor = pow(10., -0.4 * ext[j] );
         
-        for(i=0; i<cam.inumpixels;i++){
-            nsb_phepns[i]+=diffnsb_phepns[i]/iNUMWAVEBANDS + zenfactor *nsbrate_phepns[i][j];
+        for(UInt_t ui=0; ui<camgeom.GetNumPixels();ui++){
+            nsb_phepns[ui]+=diffnsb_phepns[ui]/iNUMWAVEBANDS + zenfactor *nsbrate_phepns[ui][j];
         }
     }
@@ -1327 +1327 @@
 
         k = produce_phes( inputfile,
-                          &cam,
+                          &camgeom,
                           WAVEBANDBOUND1,
                           WAVEBANDBOUND6,
@@ -1349 +1349 @@
 
           //  Fill trigger and fadc response in the trigger class from the database
-          for(i=0;i<cam.inumpixels;i++)
-            if(nsb_phepns[i]>0.0){
-              k=lons.GetResponse(nsb_phepns[i],0.1,
+          for(UInt_t ui=0;ui<camgeom.GetNumPixels();ui++)
+            if(nsb_phepns[ui]>0.0){
+              k=lons.GetResponse(nsb_phepns[ui],0.1,
                                  & nsb_trigresp[0],& nsb_fadcresp[0]);
               if(k==0){
@@ -1357 +1357 @@
                 exit(1);
               }
-              Trigger.AddNSB(i,nsb_trigresp);
-              fadc.AddSignal(i,nsb_fadcresp);
+              Trigger.AddNSB(ui,nsb_trigresp);
+              fadc.AddSignal(ui,nsb_fadcresp);
             }
         }// end if(simulateNSB && nphe2NSB<=inumphe) ...
@@ -1387 +1387 @@
         cout << "Total number of phes: " << inumphe<<" + ";
         inumphensb=0;
-        for (i=0;i<cam.inumpixels;i++){
-          inumphensb+=nsb_phepns[i]*TOTAL_TRIGGER_TIME;
+        for (UInt_t ui=0;ui<camgeom.GetNumPixels();ui++){
+          inumphensb+=nsb_phepns[ui]*TOTAL_TRIGGER_TIME;
         }
         cout<<inumphensb<<endl;
@@ -2751 +2751 @@
 
 /******** bpoint_is_in_pix() ***************************************/
-/*
-int bpoint_is_in_pix(double dx, double dy, int ipixnum, struct camera *pcam){
-  // return TRUE if point (dx, dy) is in pixel number ipixnum, else return FALSE (use camera coordinate system)
-  // the pixel is assumed to be a "closed set"
-
-  double a, b; // a = length of one of the edges of one pixel, b = half the width of one pixel
-  double c, xx, yy; // auxiliary variable
-
-  b = pcam->dpixdiameter_cm / 2. * pcam->dpixsizefactor[ipixnum];
-  a = pcam->dpixdiameter_cm / sqrt(3.) * pcam->dpixsizefactor[ipixnum];
-  c = 1./sqrt(3.);
-
-  if((ipixnum < 0)||(ipixnum >= pcam->inumpixels)){
-    fprintf(stderr, "Error in bpoint_is_in_pix: invalid pixel number %d\n", ipixnum);
-    fprintf(stderr, "Exiting.\n");
-    exit(203);
-  }
-  xx = dx - pcam->dxc[ipixnum];
-  yy = dy - pcam->dyc[ipixnum];
-
-  if(((-b <= xx) && (xx <= 0.) && ((-c * xx - a) <= yy) && (yy <= ( c * xx + a))) ||
-     ((0. <  xx) && (xx <= b ) && (( c * xx - a) <= yy) && (yy <= (-c * xx + a)))   ){
-    return(TRUE); // point is inside
-  }
-  else{
-    return(FALSE); // point is outside
-  }
-}
-*/
 
 #define sqrt13 0.577350269 // = 1./sqrt(3.)
 #define sqrt3  1.732050807 // = sqrt(3.)
 
-int bpoint_is_in_pix(double dx, double dy, struct camera *pcam)
+int bpoint_is_in_pix(double dx, double dy, MGeomCamMagic *pgeomcam)
 {
     /* return TRUE if point (dx, dy) is in pixel number ipixnum, else return FALSE (use camera coordinate system) */
@@ -2794 +2765 @@
      */
 
-    const int    num  = pcam->inumpixels;
-    const double diam = pcam->dpixdiameter_cm;
-
-    const double diam2 = diam/2;
-    const double diam3 = diam/sqrt3;
-
-    for (int i=0; i<num; i++)
+    const int    numN  = pgeomcam->GetNumPixels();
+
+    for (int i=0; i<numN; i++)
     {
-        const double size = pcam->dpixsizefactor[i];
-
-        const double b = diam2 * size;
-        const double a = diam3 * size;
-
-        const double xx = dx - pcam->dxc[i];
-        const double yy = dy - pcam->dyc[i];
-
-        if(((-b <= xx) && (xx <= 0.) && ((-sqrt13 * xx - a) <= yy) && (yy <= ( sqrt13 * xx + a))) ||
-           ((0. <  xx) && (xx <= b ) && (( sqrt13 * xx - a) <= yy) && (yy <= (-sqrt13 * xx + a)))   ){
-
-            return i; // inside i
+      MGeomPix &pixel = (*pgeomcam)[i];
+      const double b = pixel.GetR()/2;
+      const double a = pixel.GetR()/sqrt3;
+      
+      const double xx = dx - pixel.GetX();
+      const double yy = dy - pixel.GetY();
+      
+      if(((-b <= xx) && (xx <= 0.) && ((-sqrt13 * xx - a) <= yy) && (yy <= ( sqrt13 * xx + a))) ||
+         ((0. <  xx) && (xx <= b ) && (( sqrt13 * xx - a) <= yy) && (yy <= (-sqrt13 * xx + a)))   ){
+
+        return i; // inside i
         }
 
         //   return -1; // outside
     }
+
     return -1; // outside
 }
@@ -2917 +2884 @@
 
 int produce_phes( FILE *sp, // the input file
-                  struct camera *cam, // the camera layout
+                  MGeomCamMagic *camgeom, // the camera layout
                   float minwl_nm, // the minimum accepted wavelength
                   float maxwl_nm, // the maximum accepted wavelength
@@ -2930 +2897 @@
                    ){
 
-  static int i, k, ipixnum;
+  static uint i;
+  static int k, ipixnum;
   static float cx, cy, wl, qe, t;
   static MCCphoton photon;
   static float **qept;
   static char flag[SIZE_OF_FLAGS + 1];
-  static float radius;
+  static float radius_mm;
   static float SFR  , C3 , C2 , C1; 
   const double pi = 3.14159;
@@ -2941 +2909 @@
   // reset variables
 
-  for ( i=0; i<cam->inumpixels; ++i ){
+  for ( i=0; i<camgeom->GetNumPixels(); ++i ){
     
     nphe[i] = 0.0;
@@ -2949 +2917 @@
   *incph = 0;
 
-  radius = cam->dxc[cam->inumpixels-1]
-    + 1.5*cam->dpixdiameter_cm*cam->dpixsizefactor[cam->inumpixels-1];
+  radius_mm = camgeom->GetMaxRadius();
+
 
   //- - - - - - - - - - - - - - - - - - - - - - - - - 
@@ -3043 +3011 @@
     // check if photon has valid wavelength and is inside outermost camera radius
           
-    if( (wl > maxwl_nm) || (wl < minwl_nm) || (sqrt(cx*cx + cy*cy) > radius ) ){ 
+    if( (wl > maxwl_nm) || (wl < minwl_nm) || (sqrt(cx*cx + cy*cy)*10 > radius_mm ) ){ 
       continue;
             
     }
 
-    ipixnum = bpoint_is_in_pix(cx, cy, cam);
+    ipixnum = bpoint_is_in_pix(cx*10, cy*10, camgeom);
 
     // -1 = the photon is in none of the pixels
@@ -3121 +3089 @@
 
 int produce_nsbrates( char *iname, // the starfield input file name
-                      struct camera *cam, // camera layout
+                      MGeomCamMagic *camgeom, // camera layout
                       float rate_phepns[iMAXNUMPIX][iNUMWAVEBANDS] // the product of this function:
                                                // the NSB rates in phe/ns for each pixel 
                       ){
 
-  int i, j, k, ii; // counters
+  uint i, j; // counters
+  int k, ii; // counters
 
   MTrigger trigger(Trigger_gate_length, Trigger_overlaping_time, Trigger_response_ampl, Trigger_response_fwhm);
@@ -3212 +3181 @@
           // initialize the rate array
 
-          for(j = 0; j<cam->inumpixels; j++){ // loop over pixels
+          for(j = 0; j<camgeom->GetNumPixels(); j++){ // loop over pixels
             rate_phepns[j][i] = 0.;
           }
@@ -3231 +3200 @@
 
             k = produce_phes( infile,
-                              cam,
+                              camgeom,
                               wl_nm[i],
                               wl_nm[i+1],
@@ -3248 +3217 @@
             }
           
-            for(j = 0; j<cam->inumpixels; j++){ // loop over pixels
+            for(j = 0; j<camgeom->GetNumPixels(); j++){ // loop over pixels
               rate_phepns[j][i] += nphe[j]/integtime_ns/(float)iNUMNSBPRODCALLS;
             }
@@ -3409 +3378 @@
 //
 // $Log: not supported by cvs2svn $
+// Revision 1.40  2002/07/16 16:15:22  blanch
+// A first implementation of the Star field rotation has been done.
+//
 // Revision 1.39  2002/04/27 10:48:39  blanch
 // Some unused varibles have been removed.
@@ -3528 +3500 @@
 //
 // $Log: not supported by cvs2svn $
+// Revision 1.40  2002/07/16 16:15:22  blanch
+// A first implementation of the Star field rotation has been done.
+//
 // Revision 1.39  2002/04/27 10:48:39  blanch
 // Some unused varibles have been removed.