Ignore:
Timestamp:
03/16/04 12:37:24 (21 years ago)
Author:
wittek
Message:
*** empty log message ***
File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/manalysis/MSigmabar.cc

    r2798 r3522  
    100100    Int_t innerPixels[6];
    101101    Int_t outerPixels[6];
    102     Float_t innerSquaredSum[6];
    103     Float_t outerSquaredSum[6];
     102    Float_t innerSum[6];
     103    Float_t outerSum[6];
    104104
    105105    memset(innerPixels, 0, sizeof(innerPixels));
    106106    memset(outerPixels, 0, sizeof(outerPixels));
    107     memset(innerSquaredSum, 0, sizeof(innerSquaredSum));
    108     memset(outerSquaredSum, 0, sizeof(outerSquaredSum));
    109 
    110     //
    111     // sum up sigma^2 for each sector, separately for inner and outer region;
    112     // all pixels are renormalized to the area of pixel 0
     107    memset(innerSum, 0, sizeof(innerSum));
     108    memset(outerSum, 0, sizeof(outerSum));
     109
     110    //
     111    // sum up sigma/sqrt(area) for each sector,
     112    // separately for inner and outer region;
    113113    //
    114114    // consider all pixels with Cherenkov photon information
     
    149149        const MGeomPix &gpix = geom[idx];
    150150
    151         Int_t sector = (Int_t)(atan2(gpix.GetY(),gpix.GetX())*6 / (TMath::Pi()*2));
     151        Int_t sector = (Int_t)(atan2(gpix.GetY(),gpix.GetX())*6
     152                       / (TMath::Pi()*2));
    152153        if (sector<0)
    153154            sector+=6;
     
    163164        {
    164165            innerPixels[sector]++;
    165             innerSquaredSum[sector]+= sigma*sigma * ratio;
     166            innerSum[sector]+= sigma * sqrt(ratio);
    166167        }
    167168        else
    168169        {
    169170            outerPixels[sector]++;
    170             outerSquaredSum[sector]+= sigma*sigma * ratio;
     171            outerSum[sector]+= sigma * sqrt(ratio);
    171172        }
    172173    }
     
    174175    fInnerPixels   = 0;
    175176    fOuterPixels   = 0;
    176     fSigmabarInner = 0;
    177     fSigmabarOuter = 0;
     177    Double_t fSumInner = 0;
     178    Double_t fSumOuter = 0;
    178179    for (UInt_t i=0; i<6; i++)
    179180    {
    180         fSigmabarInner += innerSquaredSum[i];
     181        fSumInner      += innerSum[i];
    181182        fInnerPixels   += innerPixels[i];
    182         fSigmabarOuter += outerSquaredSum[i];
     183        fSumOuter      += outerSum[i];
    183184        fOuterPixels   += outerPixels[i];
    184185    }
    185186
    186     //
    187     // this is the sqrt of the average sigma^2;
    188     //
    189     fSigmabar=fInnerPixels+fOuterPixels<=0?0:sqrt((fSigmabarInner+fSigmabarOuter)/(fInnerPixels+fOuterPixels));
    190 
    191     if (fInnerPixels > 0) fSigmabarInner /= fInnerPixels;
    192     if (fOuterPixels > 0) fSigmabarOuter /= fOuterPixels;
    193 
    194     //
    195     // this is the sqrt of the average sigma^2
    196     // for the inner and outer pixels respectively
    197     //
    198     fSigmabarInner = sqrt( fSigmabarInner );
    199     fSigmabarOuter = sqrt( fSigmabarOuter );
     187    if (fInnerPixels > 0) fSigmabarInner = fSumInner / fInnerPixels;
     188    if (fOuterPixels > 0) fSigmabarOuter = fSumOuter / fOuterPixels;
     189
     190    //
     191    // this is the average sigma/sqrt(area) (over all pixels)
     192    //
     193    fSigmabar = (fInnerPixels+fOuterPixels)<=0 ? 0:
     194                (fSumInner+fSumOuter)/(fInnerPixels+fOuterPixels);
    200195
    201196    for (UInt_t i=0; i<6; i++)
     
    203198        const Double_t ip  = innerPixels[i];
    204199        const Double_t op  = outerPixels[i];
    205         const Double_t iss = innerSquaredSum[i];
    206         const Double_t oss = outerSquaredSum[i];
     200        const Double_t iss = innerSum[i];
     201        const Double_t oss = outerSum[i];
    207202
    208203        const Double_t sum = ip + op;
    209 
    210         fSigmabarSector[i]      = sum<=0 ? 0 : sqrt((iss+oss)/sum);
    211         fSigmabarInnerSector[i] = ip <=0 ? 0 : sqrt(iss/ip);
    212         fSigmabarOuterSector[i] = op <=0 ? 0 : sqrt(oss/op);
     204        fSigmabarInnerSector[i] = ip <=0 ? 0 :       iss/ip;
     205        fSigmabarOuterSector[i] = op <=0 ? 0 :       oss/op;
     206        fSigmabarSector[i]      = sum<=0 ? 0 : (iss+oss)/sum;
    213207    }
     208
     209    //TString opt = "";
     210    //Print(opt);
    214211
    215212  return fSigmabar;
Note: See TracChangeset for help on using the changeset viewer.