---

LayeredTube.cxx

Go to the documentation of this file.

// $Id: LayeredTube.cxx,v 1.5 2001/01/23 03:41:32 burnett Exp $
//
// This file is part of Gismo2
//

#ifdef HAVE_CONFIG_H
#include "config.h"
#endif

#include "LayeredTube.h"

#include "gismo/Detector.h"
#include "gismo/Material.h"

#include "geometry/Tube.h"
#include "geometry/Ray.h"
#include "geometry/Intersection.h"

#include "geomrep/TubeRep.h"

#include <strstream>
#include <cassert>


//                       construction

LayeredTube::LayeredTube(Medium* parent, float dz, float r)
: CompositeMedium(parent)
, m_dz(dz)
, m_r_inner(r)
, first(0), last(0)
{
}
LayeredTube::Layer::Layer(Medium* parent, const char* matname, Detector* det)
: Medium(parent, (Shape*) 0, matname, det) {}

Medium&
LayeredTube::addLayer(float thickness, const char* matname, Detector* det)
{
   Layer* newLayer = new Layer(this, matname, det);

   if( first == 0 )
   {
      first = last = newLayer;
      first->min_r = m_r_inner;
      m_r_outer = first->max_r = first->min_r+thickness;
      first->previous = 0;
      first->next = 0;
      setVolume( new Tube(m_dz, m_r_inner, m_r_outer) );
   }
   else
   {
      newLayer->min_r = last->max_r;
      m_r_outer =newLayer->max_r = newLayer->min_r+thickness;
      newLayer->previous = last;
      newLayer->next = 0;
      last->next = newLayer;
      last = newLayer;
      delete &tube();
      setVolume(new Tube(m_dz, m_r_inner, m_r_outer) );
   }

   newLayer->setVolume(new Tube(m_dz, newLayer->min_r, newLayer->max_r) );
   std::ostrstream label;
   label << "Layer " << innerMediaCount() << '\0';
   newLayer->setTitle(label.str());
   return *newLayer;
}

//                            propagation


double
LayeredTube::distanceToEnter(const Ray& r, const Medium*& m, double t)const
{
  // override this to select a Layer to enter.
  // first, do we enter the overall volume?
  double d = Medium::distanceToEnter(r,m,t);

  if( d != FLT_MAX ){

    m = select(r.position(d));
   }
  return d;

}

double
LayeredTube::distanceToLeave(const Ray& r, const Medium*& newStuff
                   , double maxDist)const
{
  newStuff = parent();
  return 0.;
}

const LayeredTube::Layer*
LayeredTube::select(const Point& x)const
{
    Vector r = x - origin();
    double rlocal = sqrt( r.mag2() - sqr(r*axis()) );
    const Layer* layer=first;

    while( (rlocal>layer->max_r) && layer->next)
    {
         layer = layer->next;
    }
    return layer;
}

const Medium*
LayeredTube::inside(const Point& x)const
{
  if( ! tube().inside(x) ) return parent();
  return select(x);
}

double
LayeredTube::Layer::distanceToLeave(const Ray& x, const Medium*& nextMed
                   , double maxStep)const
{
#ifdef _DEBUG
    const Tube* mytube=(const Tube*)(&volume());
    for( int i=0; i< mytube->surfaceCount(); i++){
        double dist = mytube->surface(i).how_near(x.position());

        assert( dist > -5e-6);
    }
#endif
    // first get distance to leave using the local tube
    double d = volume().distanceToLeave(x,  maxStep);

    if( d==FLT_MAX ) {
        std::ostrstream error_text;
        Vector r = x.position()-origin();

        error_text << "LayeredTube::Layer::distanceToLeave: point not inside??\n"
                   << "\tRay= " << x << '\n'
                   << "\tlocal radius, z = " << sqrt( r.mag2() - sqr(r*axis()) )
                   <<     ", " << r*axis() << '\n'
                   << "\tOrigin= " << origin() << "; axis= "<< axis() << '\n'
                   << "\tmin_r, max_r = " << min_r << "  " << max_r << '\n';
        error_text << '\0';
        WARNING(error_text.str());
        return d;
    }
    else if( d < maxStep) {
        // now figure out which way we are going and return the next medium
        Vector r = x.position(d)-origin();
        double local_z = r*axis();
        if( fabs(local_z)>= parentTube().length()/2.0)
            nextMed=0;
        else {
            Vector rhohat = (r - local_z*axis()).unit(); // unit vector in rho direction
            double slope = x.direction(d)*rhohat;
            nextMed = slope>0 ? next : previous;
        }
        d += 1e-7; // a little push if hit boundary
    }
    if( nextMed==0 ) {  // maybe went out the front or back
        nextMed = parent();
    }
#ifdef _DEBUG
    {
    mytube=(const Tube*)(&nextMed->volume());
    for( int i=0; i< mytube->surfaceCount(); i++){
        double dist = mytube->surface(i).how_near(x.position(d));

        assert( dist > -1e-7);
    }
    }
#endif

    return d;
}
//                              transform
GeomObject&
LayeredTube::transform(const CoordTransform& T)
{
    CompositeMedium::transform(T);
    m_origin.transform(T);
    return *this;
}


GeomObject&
LayeredTube::Layer::transform(const CoordTransform& T)
{
    Medium::transform(T);
    return *this;
}


//                              display
static const Tube* motherTube;   // file scope variable allowing layers access to enclosing box
static bool even;               // tag even/odd so can draw only every other layer
void
LayeredTube::createDetectorView(gui::DisplayRep& view)
{
    motherTube = &tube();
    even = true;
    CompositeMedium::createDetectorView(view);
}
void
LayeredTube::Layer::createDetectorView(gui::DisplayRep& view)
{
    if( even )  view.append( TubeRep( (Tube&)volume()) );
    even = ! even;
}
//                              print
void
LayeredTube::Layer::printOn(std::ostream& os)const
{
        os << " LayeredTube::Layer - title: "<<_title<<"\n"
           << "             _maxStep = " << maxStep()
           << " keCutOff = "<< kECutOff() << " \n";
//      _field->printOn(os);
        os << " Material: "
           << _material->name() << "\tr-range: "<< min_r << " to " << max_r << "\n";

 }

---