Ignore:
Timestamp:
10/09/02 19:15:28 (22 years ago)
Author:
bigongia
Message:
Version 0.6. Changed output format: added run header, changed event header,
             added ascii parameter files attached at the end of every output
             file to keep all info.
Location:
trunk/MagicSoft/Simulation/Detector/ReflectorII
Files:
10 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/Makefile

    r1431 r1535  
    1919#
    2020# $RCSfile: Makefile,v $
    21 # $Revision: 1.3 $
     21# $Revision: 1.4 $
    2222# $Author: bigongia $
    23 # $Date: 2002-07-24 14:35:45 $
     23# $Date: 2002-10-09 18:15:27 $
    2424#
    2525##################################################################
    … …  
    5555#CXXFLAGS  = -D__${SYSTEM}__ ${INCLUDES} ${OPTIM} ${DEBUG}
    5656CXXFLAGS  = -D__${SYSTEM}__ ${INCLUDES} ${OPTIM}
     57#CXXFLAGS  = -D__${SYSTEM}__ ${INCLUDES} ${DEBUG}
     58
    5759CFLAGS    = ${CXXFLAGS}
    5860FFLAGS    = ${CXXFLAGS}
    … …  
    7476        ph2cph.c \
    7577        header.c \
     78        attach.c \
    7679        reflector.c
    7780
    … …  
    9598        ph2cph.o \
    9699        header.o \
     100        attach.o \
    97101        reflector.o
    98102
    … …  
    173177# DO NOT DELETE THIS LINE -- make depend depends on it.
    174178
    175 diag.o: /usr/include/stdio.h /usr/include/features.h /usr/include/sys/cdefs.h
    176 diag.o: /usr/include/gnu/stubs.h
    177 diag.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stddef.h
    178 diag.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stdarg.h
    179 diag.o: /usr/include/bits/types.h /usr/include/libio.h
    180 diag.o: /usr/include/_G_config.h /usr/include/bits/stdio_lim.h
    181 diag.o: /usr/include/stdlib.h version.h diag.h
    182 init.o: /usr/include/stdio.h /usr/include/features.h /usr/include/sys/cdefs.h
    183 init.o: /usr/include/gnu/stubs.h
    184 init.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stddef.h
    185 init.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stdarg.h
    186 init.o: /usr/include/bits/types.h /usr/include/libio.h
    187 init.o: /usr/include/_G_config.h /usr/include/bits/stdio_lim.h
    188 init.o: /usr/include/string.h /usr/include/stdlib.h /usr/include/math.h
    189 init.o: /usr/include/bits/huge_val.h /usr/include/bits/mathdef.h
    190 init.o: /usr/include/bits/mathcalls.h version.h diag.h init.h
    191 parms.o: /usr/include/stdio.h /usr/include/features.h
    192 parms.o: /usr/include/sys/cdefs.h /usr/include/gnu/stubs.h
    193 parms.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stddef.h
    194 parms.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stdarg.h
    195 parms.o: /usr/include/bits/types.h /usr/include/libio.h
    196 parms.o: /usr/include/_G_config.h /usr/include/bits/stdio_lim.h
    197 parms.o: /usr/include/string.h /usr/include/stdlib.h /usr/include/math.h
    198 parms.o: /usr/include/bits/huge_val.h /usr/include/bits/mathdef.h
    199 parms.o: /usr/include/bits/mathcalls.h diag.h parms.h init.h
    200 geometry.o: /usr/include/stdio.h /usr/include/features.h
    201 geometry.o: /usr/include/sys/cdefs.h /usr/include/gnu/stubs.h
    202 geometry.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stddef.h
    203 geometry.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stdarg.h
    204 geometry.o: /usr/include/bits/types.h /usr/include/libio.h
    205 geometry.o: /usr/include/_G_config.h /usr/include/bits/stdio_lim.h
    206 geometry.o: /usr/include/string.h /usr/include/stdlib.h /usr/include/math.h
    207 geometry.o: /usr/include/bits/huge_val.h /usr/include/bits/mathdef.h
    208 geometry.o: /usr/include/bits/mathcalls.h diag.h geometry.h init.h
    209 atm.o: /usr/include/stdio.h /usr/include/features.h /usr/include/sys/cdefs.h
    210 atm.o: /usr/include/gnu/stubs.h
    211 atm.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stddef.h
    212 atm.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stdarg.h
    213 atm.o: /usr/include/bits/types.h /usr/include/libio.h
    214 atm.o: /usr/include/_G_config.h /usr/include/bits/stdio_lim.h
    215 atm.o: /usr/include/string.h /usr/include/math.h /usr/include/bits/huge_val.h
    216 atm.o: /usr/include/bits/mathdef.h /usr/include/bits/mathcalls.h diag.h atm.h
    217 atm.o: init.h
    218 ph2cph.o: /usr/include/stdio.h /usr/include/features.h
    219 ph2cph.o: /usr/include/sys/cdefs.h /usr/include/gnu/stubs.h
    220 ph2cph.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stddef.h
    221 ph2cph.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stdarg.h
    222 ph2cph.o: /usr/include/bits/types.h /usr/include/libio.h
    223 ph2cph.o: /usr/include/_G_config.h /usr/include/bits/stdio_lim.h
    224 ph2cph.o: /usr/include/math.h /usr/include/bits/huge_val.h
    225 ph2cph.o: /usr/include/bits/mathdef.h /usr/include/bits/mathcalls.h diag.h
    226 ph2cph.o: init.h lagrange.h
    227 header.o: /usr/include/string.h /usr/include/features.h
    228 header.o: /usr/include/sys/cdefs.h /usr/include/gnu/stubs.h
    229 header.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stddef.h
    230 header.o: header.h
    231 reflector.o: /usr/include/stdio.h /usr/include/features.h
    232 reflector.o: /usr/include/sys/cdefs.h /usr/include/gnu/stubs.h
    233 reflector.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stddef.h
    234 reflector.o: /usr/lib/gcc-lib/i386-redhat-linux/egcs-2.91.66/include/stdarg.h
    235 reflector.o: /usr/include/bits/types.h /usr/include/libio.h
    236 reflector.o: /usr/include/_G_config.h /usr/include/bits/stdio_lim.h
    237 reflector.o: /usr/include/stdlib.h /usr/include/string.h /usr/include/math.h
    238 reflector.o: /usr/include/bits/huge_val.h /usr/include/bits/mathdef.h
    239 reflector.o: /usr/include/bits/mathcalls.h version.h diag.h init.h header.h
     179diag.o: /usr/include/stdio.h /usr/include/standards.h /usr/include/sys/seek.h
     180diag.o: /usr/include/va_list.h /usr/include/sys/types.h
     181diag.o: /usr/include/mach/machine/vm_types.h /usr/include/sys/select.h
     182diag.o: /usr/include/getopt.h /usr/include/sys/limits.h
     183diag.o: /usr/include/sys/machine/machlimits.h /usr/include/sys/syslimits.h
     184diag.o: /usr/include/sys/machine/machtime.h /usr/include/sys/rt_limits.h
     185diag.o: /usr/include/stdlib.h /usr/include/stdarg.h version.h diag.h
     186init.o: /usr/include/stdio.h /usr/include/standards.h /usr/include/sys/seek.h
     187init.o: /usr/include/va_list.h /usr/include/sys/types.h
     188init.o: /usr/include/mach/machine/vm_types.h /usr/include/sys/select.h
     189init.o: /usr/include/getopt.h /usr/include/sys/limits.h
     190init.o: /usr/include/sys/machine/machlimits.h /usr/include/sys/syslimits.h
     191init.o: /usr/include/sys/machine/machtime.h /usr/include/sys/rt_limits.h
     192init.o: /usr/include/string.h /usr/include/strings.h /usr/include/stdlib.h
     193init.o: /usr/include/math.h version.h diag.h init.h
     194parms.o: /usr/include/stdio.h /usr/include/standards.h
     195parms.o: /usr/include/sys/seek.h /usr/include/va_list.h
     196parms.o: /usr/include/sys/types.h /usr/include/mach/machine/vm_types.h
     197parms.o: /usr/include/sys/select.h /usr/include/getopt.h
     198parms.o: /usr/include/sys/limits.h /usr/include/sys/machine/machlimits.h
     199parms.o: /usr/include/sys/syslimits.h /usr/include/sys/machine/machtime.h
     200parms.o: /usr/include/sys/rt_limits.h /usr/include/string.h
     201parms.o: /usr/include/strings.h /usr/include/stdlib.h /usr/include/math.h
     202parms.o: diag.h parms.h init.h
     203geometry.o: /usr/include/stdio.h /usr/include/standards.h
     204geometry.o: /usr/include/sys/seek.h /usr/include/va_list.h
     205geometry.o: /usr/include/sys/types.h /usr/include/mach/machine/vm_types.h
     206geometry.o: /usr/include/sys/select.h /usr/include/getopt.h
     207geometry.o: /usr/include/sys/limits.h /usr/include/sys/machine/machlimits.h
     208geometry.o: /usr/include/sys/syslimits.h /usr/include/sys/machine/machtime.h
     209geometry.o: /usr/include/sys/rt_limits.h /usr/include/string.h
     210geometry.o: /usr/include/strings.h /usr/include/stdlib.h /usr/include/math.h
     211geometry.o: diag.h geometry.h init.h
     212atm.o: /usr/include/stdio.h /usr/include/standards.h /usr/include/sys/seek.h
     213atm.o: /usr/include/va_list.h /usr/include/sys/types.h
     214atm.o: /usr/include/mach/machine/vm_types.h /usr/include/sys/select.h
     215atm.o: /usr/include/getopt.h /usr/include/sys/limits.h
     216atm.o: /usr/include/sys/machine/machlimits.h /usr/include/sys/syslimits.h
     217atm.o: /usr/include/sys/machine/machtime.h /usr/include/sys/rt_limits.h
     218atm.o: /usr/include/string.h /usr/include/strings.h /usr/include/math.h
     219atm.o: /usr/include/stdlib.h diag.h atm.h init.h
     220ph2cph.o: /usr/include/stdio.h /usr/include/standards.h
     221ph2cph.o: /usr/include/sys/seek.h /usr/include/va_list.h
     222ph2cph.o: /usr/include/sys/types.h /usr/include/mach/machine/vm_types.h
     223ph2cph.o: /usr/include/sys/select.h /usr/include/getopt.h
     224ph2cph.o: /usr/include/sys/limits.h /usr/include/sys/machine/machlimits.h
     225ph2cph.o: /usr/include/sys/syslimits.h /usr/include/sys/machine/machtime.h
     226ph2cph.o: /usr/include/sys/rt_limits.h /usr/include/math.h
     227ph2cph.o: /usr/include/stdlib.h diag.h init.h lagrange.h
     228header.o: /usr/include/string.h /usr/include/standards.h
     229header.o: /usr/include/sys/types.h /usr/include/mach/machine/vm_types.h
     230header.o: /usr/include/sys/select.h /usr/include/strings.h header.h
     231attach.o: /usr/include/string.h /usr/include/standards.h
     232attach.o: /usr/include/sys/types.h /usr/include/mach/machine/vm_types.h
     233attach.o: /usr/include/sys/select.h /usr/include/strings.h
     234reflector.o: /usr/include/stdio.h /usr/include/standards.h
     235reflector.o: /usr/include/sys/seek.h /usr/include/va_list.h
     236reflector.o: /usr/include/sys/types.h /usr/include/mach/machine/vm_types.h
     237reflector.o: /usr/include/sys/select.h /usr/include/getopt.h
     238reflector.o: /usr/include/sys/limits.h /usr/include/sys/machine/machlimits.h
     239reflector.o: /usr/include/sys/syslimits.h /usr/include/sys/machine/machtime.h
     240reflector.o: /usr/include/sys/rt_limits.h /usr/include/stdlib.h
     241reflector.o: /usr/include/string.h /usr/include/strings.h /usr/include/math.h
     242reflector.o: version.h diag.h init.h header.h

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/config.mk.linux

    r1432 r1535  
    1818#---------------------------------------------------------------
    1919# $RCSfile: config.mk.linux,v $
    20 # $Revision: 1.1 $
     20# $Revision: 1.2 $
    2121# $Author: bigongia $
    22 # $Date: 2002-07-24 14:43:34 $
     22# $Date: 2002-10-09 18:15:27 $
    2323##################################################################
    2424# @maintitle

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/geometry.c

    r1431 r1535  
    99extern  char line[];            /*  parsing buf. (init) */
    1010extern  char axisdev_filename[256], reflectivity_filename[256];
    11 
    12 float   mean_refl;              /*  Mirror mean reflectivity 270-610 nm.
    13                                  *  AM June 2002.
    14                                  */
    1511
    1612float   ct_Focal_mean;          /*  focal dist. (mean) (cm)     */
    … …  
    2420
    2521mirror *ct_data=NULL;           /*  ptr to mirror data          */
    26 FILE   *ct_BinaryData=NULL;     /*  binary data for mirrors     */
    27 char    ct_BinaryName[128];     /*  binary data filename        */
    2822
    2923int     nReflectivity=0;        /*  elements in refl. table     */
    … …  
    6660    FILE *datfile = fopen(datname, "r");
    6761    int current = 0;
    68 
    69     mean_refl = 0.;
    7062
    7163    if (datfile == NULL)
    … …  
    9385                             &Reflectivity[1][current]))
    9486      {
    95         // Added June 2002, AM:
    96         mean_refl += Reflectivity[1][current];
    9787        current++;
    9888        if (current >= nReflectivity) break;
    … …  
    10292
    10393    nReflectivity = current;
    104     if (current > 0)
    105       mean_refl /= (float) current;
    10694
    10795}   /*  end of ReadReflectivity  */

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/header.c

    r1431 r1535  
    22#include "header.h"
    33
    4 /*  Need to be defined the procedures to set the values in
    5     float TimeFirst;
    6     float TimeLast;     */
    7 
    8 static RflHeader rhead;  RflHeader *rheadp = &rhead;
    9 static CerHeader chead;  CerHeader *cheadp = &chead;
     4static RflEventHeader rhead;  RflEventHeader *rheadp = &rhead;
     5static CerEventHeader chead;  CerEventHeader *cheadp = &chead;
    106
    117extern float fixed_Phi, fixed_Theta;
    12 extern int   ct_NMirrors;
    13 extern float mean_refl;
     8extern int   is_Fixed_Target;
    149
    15 void TranslateHeader(RflHeader *r, CerHeader *c)
     10void TranslateHeader(RflEventHeader *r, CerEventHeader *c)
    1611{
    17   r->EvtNumber      = c->EvtNumber;
    18   r->PrimaryID      = c->PrimaryID;
    19   r->Etotal         = c->Etotal;
    20   r->Thick0         = c->Thick0;
    21   r->FirstTarget    = c->FirstTarget;
    22   r->zFirstInt      = c->zFirstInt;
    23   r->Theta          = c->Theta;
    24   r->Phi            = c->Phi;
    25   r->NumRndSeq      = c->NumRndSeq;
    26   r->RunNumber      = c->RunNumber;
    27   r->DateRun        = c->DateRun;
    28   r->Corsika_version= c->Corsika_version;
    29   r->NumObsLev      = c->NumObsLev;
    30   r->HeightLev      = c->HeightLev[0];
     12  memcpy(r,c,154*sizeof(float));
    3113
    32 
    33   r->SlopeSpec      = c->SlopeSpec;
    34   r->ELowLim        = c->ELowLim;
    35   r->EUppLim        = c->EUppLim;
    36   r->ThetaMin       = c->ThetaMin;
    37   r->ThetaMax       = c->ThetaMax;
    38   r->PhiMin         = c->PhiMin;
    39   r->PhiMax         = c->PhiMax;
    40   r->CWaveLower     = c->CWaveLower;
    41   r->CWaveUpper     = c->CWaveUpper;
    42 
    43   memcpy(r->p,         c->p,          3*sizeof(float));
    44   memcpy(r->RndData,   c->RndData,   30*sizeof(float));
    45   memcpy(r->CorePos,   c->CorePos,   40*sizeof(float));
    46 
    47   /* Next 4 variables added in June 2002, AM */
    48   r->telescopePhi      = fixed_Phi;
    49   r->telescopeTheta    = fixed_Theta;
    50   r->num_mirrors       = ct_NMirrors;
    51   r->mean_reflectivity = mean_refl;
    52 
    53   if (c->Corsika_version >= 6.)  /* Viewcone option implemented only in c6xx */
     14  if (is_Fixed_Target)  // <- Added in September 2002, AM
    5415    {
    55       if (c->viewcone_angles[0] > 0.)
    56         {
    57           printf ("ERROR: Input cer file was run using Corsika option\n");
    58           printf ("       VIEWCONE VUECON(1) VUECON(2)   with VUECON(1) > 0\n");
    59           printf ("  The reflector program only supports VUECON(1) = 0\n\n");
    60           exit(-1);
    61         }
    62       else
    63         r->ViewConeRadius = c->viewcone_angles[1]; /* degrees */
     16      /* Next 2 variables added in June 2002, AM */
     17      r->telescopePhi      = fixed_Phi;
     18      r->telescopeTheta    = fixed_Theta;
    6419    }
     20  // If no fixed target, telescope always points towards shower. AM Sep 2002
    6521  else
    66     r->ViewConeRadius = 0.;
    67 
     22    {
     23      r->telescopePhi      = cheadp->Phi;
     24      r->telescopeTheta    = cheadp->Theta;
     25    }
    6826} /*    end of TranslateHeader  */
    6927

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/header.h

    r1431 r1535  
    22#define __RFL_HEADER__
    33
    4 typedef struct                  /*  EVTH from cerfile  */
     4/*  RUNH from cerfile. See CORSIKA manual for explanations */
     5typedef struct
     6{   char RUNH[4];
     7    float RunNumber;
     8    float date;
     9    float Corsika_version;
     10    float NumObsLev;
     11    float HeightLev[10];
     12    float SlopeSpec;  /* Slope of primaries' energy spectrum */
     13    float ELowLim; 
     14    float EUppLim;    /* Limits of energy range for generation */
     15    float EGS4_flag;
     16    float NKG_flag;
     17    float Ecutoffh; 
     18    float Ecutoffm; 
     19    float Ecutoffe; 
     20    float Ecutoffg; 
     21    /* Physical constants and interaction flags (see CORSIKA manual): */
     22    float C[50];
     23    float dummy1[20]; /* not used */
     24    float CKA[40];
     25    float CETA[5];
     26    float CSTRBA[11];
     27    float dummy2[104]; /* not used */
     28    float AATM[5];
     29    float BATM[5];
     30    float CATM[5];
     31    float NFL[4];
     32}   CerRunHeader;
     33
     34
     35/*  EVTH from cerfile. See CORSIKA manual for explanations */
     36typedef struct
    537{   char        EVTH[4];
    638    float       EvtNumber;
    … …  
    71103    float       Multipl;       
    72104    float       CorePos[2][20];   
    73 
    74     float       dmmy1[14];
     105    float       SIBYLL[2];
     106    float       QGSJET[2];
     107    float       DPMJET[2];
     108    float       VENUS_cross;
     109    float       mu_mult_scat;
     110    float       NKG_range;
     111    float       EFRCTHN[2];
     112    float       WMAX[2];
     113    float       rthin_rmax;
     114
    75115    float       viewcone_angles[2]; /* (degrees) Inner and outer angles in
    76116                                      * Corsika's VIEWCONE option
    77117                                      */
    78118    float       dmmy2[119];
    79 }   CerHeader;
    80 
    81 typedef struct                  /*  EVTH from rflfile  */
    82 {   float       EvtNumber;
     119}   CerEventHeader;
     120
     121/*  EVTH from rflfile. We include first  all the information from the
     122    Corsika event header (except empty bytes) */
     123
     124typedef struct                 
     125{   char        EVTH[4];
     126    float       EvtNumber;
    83127    float       PrimaryID;
    84128    float       Etotal;   
    … …  
    97141    float       Corsika_version;
    98142
    99     float       NumObsLev;  /* Should be 1 for MAGIC simulation */
    100     float       HeightLev;  /* Observation Level */
    101 
    102     /* Changed meaning of next 9 variables. June 2002, A.Moralejo: */
    103     float       num_mirrors;
    104     float       mean_reflectivity;
     143    float       NumObsLev;
     144    float       HeightLev[10];
     145
     146    float       SlopeSpec;
     147    float       ELowLim;   
     148    float       EUppLim;   
     149
     150    float       Ecutoffh; 
     151    float       Ecutoffm; 
     152    float       Ecutoffe; 
     153    float       Ecutoffg; 
     154
     155    float       NFLAIN;
     156    float       NFLDIF;
     157    float       NFLPI0;
     158    float       NFLPIF;
     159    float       NFLCHE;
     160    float       NFRAGM;
     161 
     162    float       Bx;
     163    float       By;
     164 
     165    float       EGS4yn;
     166    float       NKGyn;
     167    float       GHEISHAyn;
     168    float       VENUSyn;
     169    float       CERENKOVyn;
     170    float       NEUTRINOyn;
     171    float       HORIZONTyn;
     172    float       COMPUTER;
     173
     174    float       ThetaMin;
     175    float       ThetaMax;
     176    float       PhiMin;
     177    float       PhiMax;
     178
     179    float       CBunchSize;
     180    float       CDetInX,CDetInY;
     181    float       CSpacInX,CSpacInY;
     182    float       CLenInX,CLenInY;
     183    float       COutput;
     184
     185    float       AngleNorthX;
     186    float       MuonInfo;
     187
     188    float       StepLength;
     189    float       CWaveLower;       
     190    float       CWaveUpper;       
     191    float       Multipl;       
     192    float       CorePos[2][20];   
     193    float       SIBYLL[2];
     194    float       QGSJET[2];
     195    float       DPMJET[2];
     196    float       VENUS_cross;
     197    float       mu_mult_scat;
     198    float       NKG_range;
     199    float       EFRCTHN[2];
     200    float       WMAX[2];
     201    float       rthin_rmax;
     202
     203    float       viewcone_angles[2];
     204    /* (degrees) Inner and outer angles in
     205     * Corsika's VIEWCONE option. This is
     206     * only possible with Corsika>6 versions. In
     207     * that case, PhiMin=PhiMax  and
     208     * ThetaMin=ThetaMax (also in this header)
     209     * indicate the axis of this cone.   
     210     */
     211
     212    /* ^^^ Up to here, the info from the CORSIKA event header. */
     213
     214    /* Telescope orientation: */
     215    float       telescopePhi;    /* rad */
     216    float       telescopeTheta;  /* rad */
     217
     218    /* Time of first and last photon: */
     219    float       TimeFirst;
     220    float       TimeLast;
     221
     222    /* 6 parameters and chi2 of the NKG fit to the longitudinal
     223     * particle distribution (see CORSIKA manual for explanation):
     224     */
    105225    float       longi_Nmax;
    106226    float       longi_t0;
    … …  
    111231    float       longi_chi2;
    112232
    113     float       SlopeSpec;
    114     float       ELowLim;   
    115     float       EUppLim;   
    116 
    117     float       ThetaMin;
    118     float       ThetaMax;
    119     float       PhiMin;
    120     float       PhiMax;
    121 
    122     float       CWaveLower;       
    123     float       CWaveUpper;       
    124     float       CorePos[2][20];   
    125     float       TimeFirst;
    126     float       TimeLast;
    127 
    128     /* AM, 23/05/2002: Changed meaning of following
    129      * three variables (which were unused before):
    130      */
    131     float       telescopePhi;    /* rad */
    132     float       telescopeTheta;  /* rad */
    133  
    134     float       ViewConeRadius; /* Degrees.
    135                                  * Radius of "view cone" when the primaries'
    136                                  * directions generated by Corsika lie within
    137                                  * a cone around a fixed direction. This is
    138                                  * only possible with Corsika>6 versions. In
    139                                  * that case, PhiMin=PhiMax  and
    140                                  * ThetaMin=ThetaMax (also in this header)
    141                                  * indicate the axis of this cone.   
    142                                  * If ViewConeRadius==0, it means that
    143                                  * the VIEWCONE option was not used.
    144                                  */
    145 
    146 
     233    /* Now the photon statistics of the event: */
    147234    float       CORSIKAPhs;     /*  Original photons written by Corsika */
    148235    float       AtmAbsPhs;      /*  Photons absorbed by the atmosphere  */
    … …  
    152239    float       OutOfChamPhs;   /*  Photons outside the camera          */
    153240    float       CPhotons;       /*  Photons reaching the camera         */
    154 }   RflHeader;
     241
     242    float       dummy[10];      /* not used */
     243}   RflEventHeader;
    155244
    156245
    … …  
    168257} Event_end;
    169258
    170 void TranslateHeader(RflHeader *r, CerHeader *c);
     259void TranslateHeader(RflEventHeader *r, CerEventHeader *c);
    171260
    172261#endif

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/init.h

    r1431 r1535  
    4848extern long  Seeds[2];                  /*  random seeds        */
    4949
    50 extern FILE *ct_BinaryData;             /*  binary data for mirrors     */
    51 extern char  ct_BinaryName[];           /*  binary data filename        */
    5250extern float ct_Focal_mean;             /*  focal dist. (mean) (cm)     */
    5351extern float ct_PSpread_mean;           /*  pt. spread fn. (mean) (cm)  */
    … …  
    102100    "    Single file version\n" \
    103101    "        - D. Bastieri & C. Bigongiari - Jan 2000\n" \
     102    "    Version 0.6\n" \
     103    "        - A. Moralejo - October 2002\n" \
    104104    "  ################################################\n\n"
    105105#define SIGN_ERROR_FTL          /*  program, version    */ \

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/parms.c

    r1431 r1535  
    1111
    1212char axisdev_filename[256], reflectivity_filename[256];
     13char geo_filename[256];
    1314
    1415/*  Prototypes  */
    … …  
    8283                if ((geofile=fopen(value_ptr, "r")) == NULL)
    8384                  FatalError(GEOF_ERROR_FTL, value_ptr);
     85                else
     86                  strcpy(geo_filename, value_ptr);
    8487                Message(GEOF_OPEN__MSG, value_ptr);
    85                 strcat(strcpy(ct_BinaryName, value_ptr), ".mirr");
    86                 ct_BinaryData = fopen(ct_BinaryName, "r");
    8788                break;
    8889              case axisdev_file:

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/ph2cph.c

    r1431 r1535  
    9494  float sx, sy;
    9595  float t1, t2;
     96  float dummy = 0.;
     97
    9698
    9799  void makeOmega(float theta, float phi);
    … …  
    255257
    256258  } else {
    257          
    258     d = (float) sqrt( b*b - 4.0*a*c );
    259          
     259
     260    /* Introduce positiveness check, AM 3/7/2002 */
     261
     262    dummy = b*b - 4.0*a*c;
     263
     264    if (dummy < 0.)           /* No intersection */
     265      return 2;
     266
     267    d = (float) sqrt(dummy);
     268
    260269    /*  two possible values for z */
    261270         
    … …  
    276285  xcut[0] = xCT[0] + rCT[0]/rCT[2]*(xcut[2]-xCT[2]);
    277286  xcut[1] = xCT[1] + rCT[1]/rCT[2]*(xcut[2]-xCT[2]);
     287     
    278288
    279289  /* CBC */
    … …  
    285295  sx = Lin2Curv( xcut[0] );
    286296  sy = Lin2Curv( xcut[1] );
    287          
     297
    288298  /* CBC */
    289299  Debug("@4 sx sy %f %f\n", sx, sy);
    … …  
    310320    if (distmirr2 < distmirr) {
    311321      i_mirror = i;
    312       distmirr = distmirr2;
     322      distmirr  = distmirr2;
    313323    }
    314324  }

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/reflector.c

    r1431 r1535  
    1313 *
    1414 * ----------------------------------------------------------------------
     15 *
     16 *  Modified: October 2002
     17 *  Author:   Abelardo Moralejo
     18 *  Version 0.6: introduced a run header and changed event header to keep
     19 *  all information from CORSIKA. In addition, now the magic.def, axisdev.dat
     20 *  and reflectivity.dat are attached at the end of the reflector output file
     21 *  to keep also all the parameters with which the program is run.
     22 *
     23 * ----------------------------------------------------------------------
    1524 */
    1625
    … …  
    3645static long writep = 0L;        /*  write ptr (on rflfile)      */
    3746
    38 extern CerHeader *cheadp;       /*  var inited in header.c      */
    39 extern RflHeader *rheadp;       /*  var inited in header.c      */
     47extern CerEventHeader *cheadp;  /*  var inited in header.c      */
     48extern RflEventHeader *rheadp;  /*  var inited in header.c      */
    4049
    4150/*  Prototypes  */
    … …  
    4352FILE *GetNextFile(char *cername);
    4453static int GetEvent(void);
    45 static int ProcessEvent(CerHeader *cheadp, FILE *cerfile, FILE *rflfile);
     54static int ProcessEvent(CerEventHeader *cheadp, FILE *cerfile, FILE *rflfile);
     55void attach(FILE *f1, FILE *f2);
     56
    4657
    4758FILE *chkf = NULL;
    … …  
    5162
    5263main(void)
    53 {   long event = 0L;            /*  event counter       */
     64{
     65    extern char axisdev_filename[256], reflectivity_filename[256];
     66    extern char geo_filename[256];
     67    FILE *dummy;
     68    long event = 0L;            /*  event counter       */
    5469
    5570    /*  Read init & geometry parms, init files and vars */
    … …  
    7186    fwrite(FLAG_END_OF_FILE, SIZE_OF_FLAGS, 1, rflfile);
    7287
    73     /*  Close file  */
     88    /* Attach parameter files used for running the program: */
     89    dummy = fopen (geo_filename,"r");
     90    attach (rflfile, dummy);
     91    fclose(dummy);
     92
     93    dummy = fopen (axisdev_filename,"r");
     94    attach (rflfile, dummy);
     95    fclose(dummy);
     96
     97    dummy = fopen (reflectivity_filename,"r");
     98    attach (rflfile, dummy);
     99    fclose(dummy);
     100
     101    /*  Close reflector output file  */
    74102    Log(RFLF_CLOSE_LOG);
    75103    fclose(rflfile);
    … …  
    82110  } /*  end of main  */
    83111
    84 static int ProcessEvent(CerHeader *cheadp, FILE *cerfile, FILE *rflfile)
     112static int ProcessEvent(CerEventHeader *cheadp, FILE *cerfile, FILE *rflfile)
    85113{   extern int absorption(float wlen, float height, float theta);
    86114    extern int ph2cph(photon *ph, cphoton *cph);
    … …  
    92120    extern void makeOmegaI(float theta, float phi);
    93121
    94     char pp[256];
    95122
    96123    /*  Various counters: phs = absphs + refphs[0..3] + cphs  */
    … …  
    172199        for (ph=0; ph<PH_IN_DATABLOCK; ph++)
    173200          {
    174             if (Photons[ph].w <= 0.) break;
     201            /* Added July 2002, AM: check integrity of photon info:
     202               Sometimes we found NaN values inside cerXXXXX file.
     203             */
     204            if (isnan(Photons[ph].w) || isnan(Photons[ph].x) ||
     205                isnan(Photons[ph].y) || isnan(Photons[ph].u) ||
     206                isnan(Photons[ph].v) || isnan(Photons[ph].t) ||
     207                isnan(Photons[ph].h))
     208              {
     209                printf("Warning: skipped one photon because its data contained Not-a-Number!\n");
     210                continue;
     211              }
     212
     213            if (! (Photons[ph].w > 0.) ) break;
    175214
    176215            CPhotons[cphs].w = Photons[ph].w;
    177216            Photons[ph].w = wlen = (float) fmod(Photons[ph].w, 1000.);
    178217
    179             /* TEMPORARY FIX, AM Nov 2001: we found that sometimes the value
     218            /* AM Nov 2001: we found that sometimes the value
    180219               stored in Photons[ph].w is not correct, and can result in a
    181220               wavelength beyond 600 nm, which makes the program crash later.
    182                Now we force wlen to its expected range:
    183                */
    184 
    185             wlen = MIN(MAX(290.,wlen),600.);
    186 
     221
     222               AM July 2002: we now simply skip the photon if the wavelength
     223               is not in the expected range, which we now take from the corsika
     224               event header (just in case we would change it in the future).
     225            */
     226
     227               if (wlen < cheadp->CWaveLower || wlen > cheadp->CWaveUpper)
     228                 {
     229                   printf("Warning: skipped one photon with strange wavelength: %f nm\n", wlen);
     230                   continue;
     231                 }
    187232
    188233            /* ADDED AM May 2002: now we take into account the telescope
    … …  
    215260                           ph2cph(&Photons[ph], &CPhotons[cphs])))
    216261              refphs[ref_type-1]++;
     262
    217263            else                /*      Photon passed   */
    218264              {
    … …  
    260306        rheadp->TimeLast  = last;
    261307
    262         /*  Update RflHeader with info on ph/cph nrs and write it  */
     308        /*  Update RflEventHeader with info on ph/cph nrs and write it  */
    263309        rheadp->CORSIKAPhs   = phs;
    264310        rheadp->AtmAbsPhs    = absphs;
    … …  
    268314        rheadp->OutOfChamPhs = refphs[3];
    269315        rheadp->CPhotons     = (long) overflow * NR_OF_CPHOTONS + cphs;
    270         fwrite(rheadp, sizeof(RflHeader), 1, rflfile);
    271 
    272 /*        for (myloop=0; myloop<sizeof(RflHeader)/4; myloop++)
     316        fwrite(rheadp, sizeof(RflEventHeader), 1, rflfile);
     317
     318/*        for (myloop=0; myloop<sizeof(RflEventHeader)/4; myloop++)
    273319 *        fprintf(chkf, "%e ", *((float *)rheadp+myloop));
    274320 *        fputc('\n', chkf);
    … …  
    284330            fseek(tmpf, 0L, SEEK_SET); /*  Start from the beginning  */
    285331            while (overflow--)
    286               { fread (CPhotons, sizeof(cphoton), NR_OF_CPHOTONS, tmpf);
    287                 fwrite(CPhotons, sizeof(cphoton), NR_OF_CPHOTONS, rflfile); }
     332              {
     333                fread (CPhotons, sizeof(cphoton), NR_OF_CPHOTONS, tmpf);
     334                fwrite(CPhotons, sizeof(cphoton), NR_OF_CPHOTONS, rflfile);
     335              }
    288336
    289337            /*  Reload data in CPhotons    */
    … …  
    291339
    292340            /*  Close (and delete) temp file  */
    293             fclose(tmpf);  }
     341            fclose(tmpf); 
     342          }
    294343
    295344        /*  Write (remaining) cphotons  */
     345
    296346        fwrite(CPhotons, sizeof(cphoton), cphs, rflfile);
    297347
    … …  
    308358
    309359static int GetEvent(void)
    310 {   int found = FALSE,          /*  event found         */
    311       isWrong = FALSE;          /*  cerfile is wrong    */
    312       static int newFile = TRUE; /*  if TRUE, check if cerfile is valid */
    313 
    314       do
    315         { /* In case the just-opened file is a valid cerfile,
    316              starting with a RUNH, loop until a valid event is found:
    317              id est with: first_Event <= EvtNumber <= last_Event
    318              and             low_Ecut <=  Etotal   <= high_Ecut
    319              If the search was successful, "found" is set to TRUE.
    320              If there are reading errors, "isWrong" is set to TRUE. */
     360{   
     361  int found = FALSE,          /*  event found           */
     362  isWrong = FALSE;            /*  cerfile is wrong      */
     363  static int newFile = TRUE;  /*  if TRUE, check if cerfile is valid    */
     364
     365  CerRunHeader cRunHead;
     366 
     367  do
     368    {
     369      /* In case the just-opened file is a valid cerfile,
     370         starting with a RUNH, loop until a valid event is found:
     371         id est with: first_Event <= EvtNumber <= last_Event
     372         and             low_Ecut <=  Etotal   <= high_Ecut
     373         If the search was successful, "found" is set to TRUE.
     374         If there are reading errors, "isWrong" is set to TRUE. */
    321375
    322376        if (newFile
    323          && (1 != fread(cheadp, sizeof(CerHeader), 1, cerfile)
     377         && (1 != fread(cheadp, sizeof(CerEventHeader), 1, cerfile)
    324378             || strncmp(cheadp->EVTH, "RUNH", 4)))
    325379        {   isWrong = TRUE;   }
    326380
    327381        else do
    328         {   /*  Push fileptr in case it is a "EVTH"  */
    329             readp = ftell(cerfile);
    330 
    331             /*  Error: exit loop  */
    332             if (1 != fread(cheadp, sizeof(CerHeader), 1, cerfile))
     382        {
     383          if (newFile)
     384            {
     385              /* Write Reflector "run header" (one per cer file!): */
     386              memcpy(&cRunHead, cheadp, sizeof(CerEventHeader));
     387              fwrite(&cRunHead, sizeof(CerRunHeader), 1, rflfile);
     388            }
     389
     390          /*  Push fileptr in case it is a "EVTH"  */
     391          readp = ftell(cerfile);
     392
     393          /*  Error: exit loop  */
     394          if (1 != fread(cheadp, sizeof(CerEventHeader), 1, cerfile))
    333395            {   isWrong = TRUE;    break;   }
    334396
    335             /*  Ok: set found at TRUE and exit loop  */
    336             if (strncmp(cheadp->EVTH, "EVTH", 4) == 0
    337              && first_Event <= (long)cheadp->EvtNumber
    338              &&                (long)cheadp->EvtNumber <= last_Event
    339              &&    low_Ecut <= cheadp->Etotal
    340              &&                cheadp->Etotal <= high_Ecut)
     397          /*  Ok: set found at TRUE and exit loop  */
     398          if (strncmp(cheadp->EVTH, "EVTH", 4) == 0
     399              && first_Event <= (long)cheadp->EvtNumber
     400              &&                       (long)cheadp->EvtNumber <= last_Event
     401              &&    low_Ecut <= cheadp->Etotal
     402              &&                       cheadp->Etotal <= high_Ecut)
    341403            {   found = TRUE;   }
    342404
    343             /*  File is finished: exit loop  */
    344             else if (strncmp(cheadp->EVTH, "RUNE", 4) == 0)
    345                 break;
     405          /*  File is finished: exit loop  */
     406          else if (strncmp(cheadp->EVTH, "RUNE", 4) == 0)
     407            break;
    346408
    347409        } while(found == FALSE);

  • trunk/MagicSoft/Simulation/Detector/ReflectorII/version.h

    r1431 r1535  
    66
    77#define PROGRAM reflector
    8 #define VERSION 0.5
     8#define VERSION 0.6
    99
    1010#endif
Note: See TracChangeset for help on using the changeset viewer.