source: trunk/MagicSoft/Mars/mgeom/MGeomCamDwarf.cc@ 9235

/* ======================================================================== *\
!
! *
! * This file is part of MARS, the MAGIC Analysis and Reconstruction
! * Software. It is distributed to you in the hope that it can be a useful
! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
! * It is distributed WITHOUT ANY WARRANTY.
! *
! * Permission to use, copy, modify and distribute this software and its
! * documentation for any purpose is hereby granted without fee,
! * provided that the above copyright notice appear in all copies and
! * that both that copyright notice and this permission notice appear
! * in supporting documentation. It is provided "as is" without express
! * or implied warranty.
! *
!
!
!   Author(s): Thomas Bretz 03/2007 <mailto:tbretz@astro.uni-wuerzburg.de>
!   Author(s): Michael Backes 03/2007 <mailto:michael.backes@udo.edu>
!
!   Copyright: MAGIC Software Development, 2000-2008
!
!
\* ======================================================================== */

/////////////////////////////////////////////////////////////////////////////
//
// MGeomCamDwarf
//
// This class stores the geometry information of the Dwarf camera.
//    MGeomCamDwarf cam;        // Creates the 313 pixel dwarf camera
//
// It can also be used to create a hexagonal camera with identical sized
// pixels and n rings (while the central pixel is counted as ring 0).
//    MGeomCamDwarf cam(9, 21); // Creates the CT3 camera
//
// Or it can be used to create a roundish camera, similar to a
// hexagonal camera, but the edges filled with additional pixels
// inside a circle.
//    MGeomCamDwarf cam(
//
////////////////////////////////////////////////////////////////////////////
#include "MGeomCamDwarf.h"

#include <iostream>

#include <TMath.h>
#include <TArrayI.h>

#include "MGeomPix.h"

ClassImp(MGeomCamDwarf);

using namespace std;

// --------------------------------------------------------------------------
//
//  Dwarf camera has 313 pixels. For geometry and Next Neighbor info see
//  CreateCam and CreateNN
//
MGeomCamDwarf::MGeomCamDwarf(const char *name)
    : MGeomCam(CalcNumPix(9.5), 4.57, name, "Geometry information of Dwarf Camera")
{
    CreateCam(21, 9.5);
    CreateNN();
    InitGeometry();
}

// --------------------------------------------------------------------------
//
//  Use this to create a camera with a roundish shape and a radius rad in
// millimeter containing the pixel centers. The pixel will have a diameter
// diameter in millimeters, and a distance dist in meters.
//
MGeomCamDwarf::MGeomCamDwarf(Double_t rad, Double_t diameter, Double_t dist, const char *name)
    : MGeomCam(CalcNumPix(diameter<=0 ? rad : rad/diameter), dist, name, "Geometry information for a roundish camera")
{
    CreateCam(diameter, diameter<=0 ? rad : rad/diameter);
    CreateNN();
    InitGeometry();
}

// --------------------------------------------------------------------------
//
//  Use this to create a camera with a hexagonal shape and rings rings.
// The first ring around the central pixel is 1. The pixel will have a
// diameter diameter in millimeters, and a distance dist in meters.
//
MGeomCamDwarf::MGeomCamDwarf(Int_t rings, Double_t diameter, Double_t dist, const char *name)
    : MGeomCam(CalcNumPix(rings), dist, name, "Geometry information for a hexagonal camera")
{
    CreateCam(diameter, rings);
    CreateNN();
    InitGeometry();
}

// --------------------------------------------------------------------------
//
//  Create a clone of this container. This is very easy, because we
//  simply have to create a new object of the same type.
//
TObject *MGeomCamDwarf::Clone(const char *newname) const
{
    MGeomCam *cam = new MGeomCam(GetNumPixels(), GetCameraDist());
    for (UInt_t i=0; i<GetNumPixels(); i++)
        (*this)[i].Copy((*cam)[i]);

    cam->InitGeometry();
    return cam;
}

// --------------------------------------------------------------------------
//
// Check if the photon which is flying along the trajectory u has passed
// (or will pass) the frame of the camera (and consequently get
// absorbed). The position p and direction u must be in the
// telescope coordinate frame, which is z parallel to the focal plane,
// x to the right and y upwards, looking from the mirror towards the camera.
//
// The units are cm.
//
Bool_t MGeomCamDwarf::HitFrame(MQuaternion p, const MQuaternion &u) const
{
    // z is defined from the mirror (0) to the camera (z>0).
    // Thus we just propagate to the focal plane (z=fDist)
    //p -= 1700./u.Z()*u;
    p.PropagateZ(u, GetCameraDist()*100); // m->cm

    // Add 10% to the max radius and convert from mm to cm
    return p.R()<GetMaxRadius()*0.11;//TMath::Abs(p.X())<65 && TMath::Abs(p.Y())<65;
}

// --------------------------------------------------------------------------
//
// Check if the position given in the focal plane (so z can be ignored)
// is a position which might hit the detector. It is meant to be a rough
// and fast estimate not a precise calculation. All positions dicarded
// must not hit the detector. All positions accepted might still miss
// the detector.
//
// The units are cm.
//
Bool_t MGeomCamDwarf::HitDetector(const MQuaternion &v) const
{
    // Add 10% to the max radius and convert from mm to cm
    return v.R()<GetMaxRadius()*0.11;//60.2;
}

// --------------------------------------------------------------------------
//
// Calculate in the direction 0-5 (kind of sector) in the ring-th ring
// the x and y coordinate of the i-th pixel. The unitx are unity,
// distance to (0,0) is retruned.
//
Double_t MGeomCamDwarf::CalcXY(Int_t dir, Int_t ring, Int_t i, Double_t &x, Double_t &y)
{
    const Double_t kSqrt32 = MGeomPix::gsTan60/2;

    switch (dir)
    {
    case kDirCenter: // Center
        x = 0;
        y = 0;
        break;

    case kDirNE: // Direction North East
        x = ring-i*0.5;
        y = i*kSqrt32;
        break;

    case kDirN: // Direction North
        x = ring*0.5-i;
        y = ring*kSqrt32;
        break;

    case kDirNW: // Direction North West
        x = -(ring+i)*0.5;
        y = (ring-i)*kSqrt32;
        break;

    case kDirSW: // Direction South West
         x = 0.5*i-ring;
         y = -i*kSqrt32;
         break;

    case kDirS: // Direction South
         x = i-ring*0.5;
         y = -ring*kSqrt32;
         break;

    case kDirSE: // Direction South East
        x = (ring+i)*0.5;
        y = (-ring+i)*kSqrt32;
        break;
    }
    return TMath::Hypot(x, y);
}

// --------------------------------------------------------------------------
//
// Calculate in the direction 0-5 (kind of sector) in the ring-th ring
// the x and y coordinate of the i-th pixel. The units are unity,
// distance to (0,0) is retruned.
//
Int_t MGeomCamDwarf::CalcNumPix(Int_t rings)
{
    //
    //  add the first pixel to the list
    //
    Int_t cnt = 1;

    for (Int_t ring=0; ring<rings; ring++)
        cnt += 6*(ring+1);

    return cnt;
}

// --------------------------------------------------------------------------
//
// Calculate in the direction 0-5 (kind of sector) in the ring-th ring
// the x and y coordinate of the i-th pixel. The unitx are unity,
// distance to (0,0) is retruned.
//
// Due to possible rounding errors we need to use exactly the same
// algorithm as for creating the pixels!
//
Int_t MGeomCamDwarf::CalcNumPix(Double_t rad)
{
    //
    //  add the first pixel to the list
    //
    Int_t cnt  = 1;
    Int_t ring = 1;
    while (1)
    {
        Int_t n = 0;

        //
        // calc. coords for this ring counting from the
        // starting number to the ending number
        //
        for (Int_t dir=kDirNE; dir<=kDirSE; dir++)
        {
            for (int i=0; i<ring; i++)
            {
                Double_t x, y;
                if (CalcXY(dir, ring, i, x, y)<rad)
                    n++;
            }
        }

        if (n==0)
            return cnt;

        ring++;
        cnt += n;
    }

    return cnt;
}


// --------------------------------------------------------------------------
//
//  This fills the geometry information for a hexagonal camera
//
void MGeomCamDwarf::CreateCam(Double_t diameter, Int_t rings)
{
    //    units for diameter are mm

    //
    //  add the first pixel to the list
    //
    (*this)[0].Set(0, 0, diameter);

    Int_t cnt  = 1;

    for (int ring=1; ring<=rings; ring++)
    {
        for (Int_t dir=kDirNE; dir<=kDirSE; dir++)
        {
            for (int i=0; i<ring; i++)
            {
                Double_t x, y;
                CalcXY(dir, ring, i, x, y);
                (*this)[cnt++].Set(x*diameter, y*diameter, diameter);
            }
        }
    }
}

// --------------------------------------------------------------------------
//
//  This fills the geometry information for a roundish camera
//
void MGeomCamDwarf::CreateCam(Double_t diameter, Double_t rad)
{
    //    units for diameter are mm

    //
    //  add the first pixel to the list
    //
    (*this)[0].Set(0, 0, diameter);

    Int_t cnt  = 1;
    Int_t ring = 1;

    while (1)
    {
        Int_t n = 0;

        for (Int_t dir=kDirNE; dir<=kDirSE; dir++)
        {
            for (int i=0; i<ring; i++)
            {
                Double_t x, y;
                if (CalcXY(dir, ring, i, x, y)<rad)
                    (*this)[cnt+n++].Set(x*diameter, y*diameter, diameter);
            }
        }

        if (n==0)
            return;

        ring++;
        cnt += n;
    }
}

// --------------------------------------------------------------------------
//
//  This fills the next neighbor information from a table into the pixel
//  objects.
//
void MGeomCamDwarf::CreateNN()
{
    TArrayI nn(6);

    for (UInt_t i=0; i<GetNumPixels(); i++)
    {
        MGeomPix &pix = (*this)[i];

        Int_t k = 0;
        nn.Reset(-1);

        for (UInt_t j=0; j<GetNumPixels(); j++)
        {
            if (i==j)
                continue;

            if (pix.GetDist((*this)[j])<pix.GetD()*1.5)
                nn[k++] = j;
        }

        pix.SetNeighbors(nn[0], nn[1], nn[2], nn[3], nn[4], nn[5]);
    }
}