MdcGeomSvc-00-01-42/src/MdcGeoWire.cxx

Go to the documentation of this file.

 
 
 
#include "MdcGeomSvc/MdcGeoWire.h"
#include "MdcGeomSvc/MdcGeomSvc.h"
#include <math.h>
 
// according to the mathematicial formula:
// Sg = L^2*g*rhol/8*T
// rhol is the linear density 
// rhol is fixed at 9.47E-3 g/m
const double MdcGeoWire::Sag(void) const{
 
    if(MdcGeomSvc::getSagFlag()){
    // unit of wire length is mm
    const double L = (fForward - fBackward).mag(); 
    // sag unit is mm
    return (L*L*1.184E-6/fTension);
    }
    // do not consider wire sagita
    else return 0.;
}
 
 
// calculate the sag at length z
// attention, the z is not the z-coordinate of the point
// but the length from this point to the backward point 
const double MdcGeoWire::Sagz(const double z) const{
  return z*0.; //FIXME
}
 
 
void MdcGeoWire::AddWireNeighbor(int id,int type){
    fWireNeighborsID.push_back(id);
    fNeighborType.push_back(type);
}
 
void MdcGeoWire::AddWireNeighbor(MdcGeoWire* awire,int type){
    fWireNeighbors.push_back(awire);
    fNeighborType.push_back(type);
}
 
void MdcGeoWire::AddWireNeighbor(MdcGeoWire* awire,int id,int type){
     fWireNeighbors.push_back(awire);
    fNeighborType.push_back(type);
    fWireNeighborsID.push_back(id);
}
 
double MdcGeoWire::GetX(double z) const{
    return (fForward.x()-fBackward.x())*(z-fBackward.z())/(fForward.z()-fBackward.z())+fBackward.x();
}
 
double MdcGeoWire::GetY(double z) const{
     return (fForward.y()-fBackward.y())*(z-fBackward.z())/(fForward.z()-fBackward.z())+fBackward.y();
}
 
vector<int> MdcGeoWire::GetNeighborIDType3(void) const{
    std::vector<int> NeighborType3;
    for(unsigned i=0;i<fWireNeighborsID.size();i++){
        if(fNeighborType.at(i)==3){
            NeighborType3.push_back(fWireNeighborsID.at(i));
        }
    }
    return NeighborType3;
}
 
vector<int> MdcGeoWire::GetNeighborIDType4(void) const{
    std::vector<int> NeighborType4;
    for(unsigned i=0;i<fWireNeighborsID.size();i++){
        if(fNeighborType.at(i)==4){
            NeighborType4.push_back(fWireNeighborsID.at(i));
        }
    }
    return NeighborType4;
}
 
vector<MdcGeoWire*> MdcGeoWire::GetNeighborType3(void) const{
    std::vector<MdcGeoWire*> NeighborType3;
    for(unsigned i=0;i<fWireNeighbors.size();i++){
        if(fNeighborType.at(i)==3){
            NeighborType3.push_back(fWireNeighbors.at(i));
        }
    }
    return NeighborType3;
}
 
vector<MdcGeoWire*> MdcGeoWire::GetNeighborType4(void) const{
    std::vector<MdcGeoWire*> NeighborType4;
    for(unsigned i=0;i<fWireNeighbors.size();i++){
        if(fNeighborType.at(i)==4){
            NeighborType4.push_back(fWireNeighbors.at(i));
        }
    }
    return NeighborType4;
}
 
 
vector<MdcGeoWire*> MdcGeoWire::GetNeighboratZType3(double z) const{
    std::vector<MdcGeoWire*> NeighborZType3;
    const double phi = MdcGeoWire::WirePhi(z);
    double phi2,dPhi1,dPhi2,phiMax,phiMin;
    std::vector<MdcGeoWire*> NeighborType3 = MdcGeoWire::GetNeighborType3();
    dPhi1=M_PI/double(fLyr->NCell());
    if(NeighborType3.size()>0) {
        dPhi2=M_PI/double(NeighborType3.at(0)->Lyr()->NCell());
        phiMax=MdcGeoWire::CalculatePhi(phi+dPhi1+dPhi2);
        phiMin=MdcGeoWire::CalculatePhi(phi-dPhi1-dPhi2);
        if(fabs(z)<NeighborType3.at(0)->Lyr()->Length()/2.){
            for(unsigned i=0;i<NeighborType3.size();i++){
                phi2 =NeighborType3.at(i)->WirePhi(z);
                if(phi2>phiMin&&phi2<phiMax) NeighborZType3.push_back(NeighborType3.at(i));
                else if(phiMin>phiMax&&(phi2>phiMin||phi2<phiMax)) NeighborZType3.push_back(NeighborType3.at(i));
            }
        }
    }
    return NeighborZType3;
}
 
vector<MdcGeoWire*> MdcGeoWire::GetNeighboratZType4(double z) const{
    std::vector<MdcGeoWire*> NeighborZType4;
    const double phi = MdcGeoWire::WirePhi(z);
    double phi2,dPhi1,dPhi2,phiMax,phiMin;
    std::vector<MdcGeoWire*> NeighborType4 = MdcGeoWire::GetNeighborType4();
    dPhi1=M_PI/double(fLyr->NCell());
    if(NeighborType4.size()>0) {
        dPhi2=M_PI/double(NeighborType4.at(0)->Lyr()->NCell());
        phiMax=MdcGeoWire::CalculatePhi(phi+dPhi1+dPhi2);
        phiMin=MdcGeoWire::CalculatePhi(phi-dPhi1-dPhi2);
        
        for(unsigned i=0;i<NeighborType4.size();i++){
            phi2 =NeighborType4.at(i)->WirePhi(z);
            if(phi2>phiMin&&phi2<phiMax) NeighborZType4.push_back(NeighborType4.at(i));
            else if(phiMin>phiMax&&(phi2>phiMin||phi2<phiMax)) NeighborZType4.push_back(NeighborType4.at(i));
        }
 
    }
    return NeighborZType4;
}
 
void MdcGeoWire::ClearNeighbors(){
  fWireNeighborsID.clear();
  fWireNeighbors.clear();
  fNeighborType.clear();
}
 
double MdcGeoWire::WirePhi(double z) const{
    double x,y;
    double phi;
    x = MdcGeoWire::GetX(z);
    y = MdcGeoWire::GetY(z);
    if(y>=0){
            phi = acos(x/sqrt(pow(x,2.)+pow(y,2.)));
    }
    else {
            phi =2*M_PI-acos(x/sqrt(pow(x,2.)+pow(y,2.)));
    }
    return phi;
}
 
double MdcGeoWire::CalculatePhi(double phi) const{
     if(0<phi&&phi<2*M_PI){
         return phi;
     }
     else if(phi>2*M_PI){
         return phi-2*M_PI;
     }
     else {
         return phi+2*M_PI;
     }
 }