[u/mrichter/AliRoot.git] / ITS / AliITSBaseGeometry.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/
/*
  $Id:
 */

#include <TObject.h>
#include <TObjArray.h>

#include <TGeoMaterial.h>
#include <TGeoMatrix.h>
#include <TGeoMedium.h>
#include <TGeoPgon.h>
#include "AliITSBaseGeometry.h"


ClassImp(AliITSMixture)

AliITSMixture::AliITSMixture(const char *name,Int_t N,Double_t *w,TObjArray *m,
			     Double_t rho,Double_t radlen,Double_t intleng)
    :TGeoMixture(name,1,rho){
    // Defines a new mixture from a number of Mixtures, and put the
    // resulting mixture into this object. This will compute avarage
    // isotopic value between different elements.
    // Inputs:
    //   Int_t    N   The number of mixtures in te TObjArray
    //   Double_t *w  The array of weights of each mixture
    //   TObjArray *m The array of AliITSMixture (TGeoMixture)s 
    //                to be mixed.
    // Output:
    //   none.
    // Return:
    //   none.
    Int_t i,z=0,j,Nel;
    Double_t tw,*nw,wel[110],Ael[110],el[110];
    TGeoMixture *mix;

    if(N>m->GetEntries()){ // Error not enough mixtures defined
	Error("Mixing","There are more weight defined than mixtures");
	return;
    } // end if
    // First normilize the weights just in case.
    tw = 0.0;
    for(i=0;i<N;i++) if(w[i]>0.0) tw += w[i];
    nw = new Double_t[N];
    for(i=0;i<N;i++) {if(w[i]>0.0) nw[i] = w[i]/tw;else nw[i] = 0.0;}
    //
    Nel=0;
    for(i=0;i<110;i++) {el[i] = wel[i] = Ael[i] = 0.0;}
    for(i=0;i<N;i++)if(w[i]>0.0) {
	mix = (TGeoMixture*) (m->At(i));
	for(j=0;j<mix->GetNelements();j++) {
	    z = (Int_t) ((mix->GetZmixt())[j]);
	    wel[z] += nw[i]*((mix->GetWmixt())[j]);
	    el[z]  += wel[z]*((mix->GetZmixt())[j]);
	    Ael[z] += wel[z]*((mix->GetAmixt())[j]);
	} // end for j
    } // end for i
    tw = 0.0;
    for(i=1;i<110;i++) if(wel[i]>0.0){
	Nel++;
	tw += wel[i];
    } // end for
    if(tw<=0.0) {  // Error no elements defined.
	Error("Mixing","Total weight of this mixture is zero");
	delete[] nw;
	return;
    } // end if
    // setup TGeoMixture data members.
    fNelements = Nel;
    if(fZmixture!=0) delete[] fZmixture;
    if(fAmixture!=0) delete[] fAmixture;
    if(fWeights!=0)  delete[] fWeights;
    fZmixture = new Double_t[Nel];
    fAmixture = new Double_t[Nel];
    fWeights  = new Double_t[Nel];
    if(rho>0.) fDensity = rho;
    else { // try to compute density form mixture.
	rho = 0.0;
	for(i=0;i<N;i++) if(nw[i]>0.0) {
	    mix = (TGeoMixture*) (m->At(i));
	    rho += nw[i]*(mix->GetDensity());
	} // end for i
	fDensity = rho;
    } // end if
    if(radlen>0.) fRadLen = radlen;
    else { // try to compute radiation form mixture.
	// From "Review of Particle Physics" Particle Data Group Section
	// 26.4.1 equation 26.21 (2002).
	radlen = 0.0;
	for(i=0;i<N;i++) if(nw[i]>0.0) {
	    mix = (TGeoMixture*) (m->At(i));
	    if(mix->GetRadLen()>0.0) rho += 1.0/(nw[i]*(mix->GetRadLen()));
	} // end for i
	fRadLen = 1.0/radlen;
    } // end if
    if(intleng>0.) fIntLen = intleng;
    else { // try to compute interaction form mixture.
	intleng = 0.0;
	for(i=0;i<N;i++) if(nw[i]>0.0) {
	    mix = (TGeoMixture*) (m->At(i));
	    if(mix->GetIntLen()>0.0) intleng += 1.0/(nw[i]*(mix->GetIntLen()));
	} // end for i
	fIntLen = 1.0/intleng;
    } // end if
    j = 0;
    for(z=1;z<110;z++){
	wel[z] /= tw;
	el[z]  /= tw;
	Ael[z] /= tw;
	if(wel[z]>0.0) this->DefineElement(j++,Ael[z],el[z],wel[z]);
    } // end for i
    delete[] nw;
}
//======================================================================
ClassImp(AliITSGeoCable)
;
AliITSGeoCable::AliITSGeoCable(){
    //

    fRmin = fRmax = 0.0;
    fNs.SetXYZ(0.0,0.0,0.0);
    fNe.SetXYZ(0.0,0.0,0.0);
    fTubes    = 0;
    fTranRot  = 0;
}
//----------------------------------------------------------------------
AliITSGeoCable::AliITSGeoCable(const char *name,const TObjArray *vect,
                               const Double_t Rmin,const Double_t Rmax,
                               const TVector3 ns,const TVector3 ne){
    //
    // Inputs:
    //    char      *name   Name of this compound object
    //    TObjArray *vect   Array of TVector3's of points representing the
    //                      path of the cable
    //    TVector3  ns=0    Normal vector representing the angle of the 
    //                      starting surface, default perpendicular
    //    TVector3  ne=0    Normal vector representing the angle of the 
    //                      ending surface, default perpendicular
    // Outputs:
    //    none.
    // Return:
    //    A fully initilized and created AliITSGeoCable class.
    Char_t nam[500];
    Int_t i,n;
    Double_t s,th,ph;
    TVector3 x0,x1,x2,d,t,n0,n1;

    fRmin = Rmin;
    fRmax = Rmax;
    fRmin = fRmax = 0.0;
    fNs.SetXYZ(0.0,0.0,-1.0);
    fNe.SetXYZ(0.0,0.0,1.0);
    n = vect->GetEntries();
    fTubes = new TObjArray(n-1);
    fTranRot  = new TObjArray(n-1);
    fNs = ns*(1./ns.Mag());
    fNe = ne*(1./ne.Mag());
    //
    x0 = *((TVector3 *)(vect->At(0)));
    n0 = ns;
    for(i=1;i<n;i++){
        x1 = *((TVector3 *)(vect->At(i)));
        d  = x1 - x0;
        if(i<n-1) {
            x2 = *((TVector3 *)(vect->At(i+1)));
            n1 = d + (x2-x1);
            n1 *= 1./n1.Mag();
        }else{
            n1 = fNe;
        } // end if
        t  = 0.5*(x1 + x0);
        th =  TMath::ATan2(TMath::Sqrt(2.*d.Mag2()-d.z()-2.*d.Mag()*d.y()-
                                        2.*d.Mag()*d.x()),
                           TMath::Sqrt(d.z()*d.z()-2.*d.Mag()*d.z()+d.Mag2()));
        th *= TMath::RadToDeg();
        ph =  TMath::ATan2(d.y()-d.Mag(),d.x()-d.Mag());
        ph *= TMath::RadToDeg();
        sprintf(nam,"%sCombiTrans%dCable",name,i-1);
        fTranRot->AddAt(new TGeoCombiTrans(nam,t.x(),t.y(),t.z(),
					 new TGeoRotation("",ph,th,0.0)),i-1);
        s  = d.Mag();
        sprintf(nam,"%sPart%dCable",name,i-1);
        fTubes->AddAt( new TGeoCtub(nam,fRmin,fRmax,0.5*s,0.0,360.0,n0.x(),
				    n0.y(),n0.z(),n1.x(),n1.y(),n1.z()),i-1);
        n0 = -n1;
        x0 = x1;
    } // end for i
}
//----------------------------------------------------------------------
AliITSGeoCable::~AliITSGeoCable(){
    //
    // Inputs:
    //    none.
    // Outputs:
    //    none.
    // Return:
    //    none.
    Int_t i;

    if(fTubes){
        for(i=0;i<fTubes->GetEntries();i++) 
            delete (TGeoCtub*)(fTubes->At(i));
        delete fTubes;
    } // end if
    fTubes = 0;
    if(fTranRot){
        for(i=0;i<fTranRot->GetEntries();i++) 
            delete (TGeoCombiTrans*)(fTranRot->At(i));
        delete fTranRot;
    } // end if
    fTranRot  = 0;
}
//----------------------------------------------------------------------
Commit	Line	Data
	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15	/*
	16	$Id:
	17	*/
	18
	19	#include <TObject.h>
	20	#include <TObjArray.h>
	21
	22	#include <TGeoMaterial.h>
	23	#include <TGeoMatrix.h>
	24	#include <TGeoMedium.h>
	25	#include <TGeoPgon.h>
	26	#include "AliITSBaseGeometry.h"
	27
	28
	29	ClassImp(AliITSMixture)
	30
	31	AliITSMixture::AliITSMixture(const char name,Int_t N,Double_t w,TObjArray *m,
	32	Double_t rho,Double_t radlen,Double_t intleng)
	33	:TGeoMixture(name,1,rho){
	34	// Defines a new mixture from a number of Mixtures, and put the
	35	// resulting mixture into this object. This will compute avarage
	36	// isotopic value between different elements.
	37	// Inputs:
	38	// Int_t N The number of mixtures in te TObjArray
	39	// Double_t *w The array of weights of each mixture
	40	// TObjArray *m The array of AliITSMixture (TGeoMixture)s
	41	// to be mixed.
	42	// Output:
	43	// none.
	44	// Return:
	45	// none.
	46	Int_t i,z=0,j,Nel;
	47	Double_t tw,*nw,wel[110],Ael[110],el[110];
	48	TGeoMixture *mix;
	49
	50	if(N>m->GetEntries()){ // Error not enough mixtures defined
	51	Error("Mixing","There are more weight defined than mixtures");
	52	return;
	53	} // end if
	54	// First normilize the weights just in case.
	55	tw = 0.0;
	56	for(i=0;i<N;i++) if(w[i]>0.0) tw += w[i];
	57	nw = new Double_t[N];
	58	for(i=0;i<N;i++) {if(w[i]>0.0) nw[i] = w[i]/tw;else nw[i] = 0.0;}
	59	//
	60	Nel=0;
	61	for(i=0;i<110;i++) {el[i] = wel[i] = Ael[i] = 0.0;}
	62	for(i=0;i<N;i++)if(w[i]>0.0) {
	63	mix = (TGeoMixture*) (m->At(i));
	64	for(j=0;j<mix->GetNelements();j++) {
	65	z = (Int_t) ((mix->GetZmixt())[j]);
	66	wel[z] += nw[i]*((mix->GetWmixt())[j]);
	67	el[z] += wel[z]*((mix->GetZmixt())[j]);
	68	Ael[z] += wel[z]*((mix->GetAmixt())[j]);
	69	} // end for j
	70	} // end for i
	71	tw = 0.0;
	72	for(i=1;i<110;i++) if(wel[i]>0.0){
	73	Nel++;
	74	tw += wel[i];
	75	} // end for
	76	if(tw<=0.0) { // Error no elements defined.
	77	Error("Mixing","Total weight of this mixture is zero");
	78	delete[] nw;
	79	return;
	80	} // end if
	81	// setup TGeoMixture data members.
	82	fNelements = Nel;
	83	if(fZmixture!=0) delete[] fZmixture;
	84	if(fAmixture!=0) delete[] fAmixture;
	85	if(fWeights!=0) delete[] fWeights;
	86	fZmixture = new Double_t[Nel];
	87	fAmixture = new Double_t[Nel];
	88	fWeights = new Double_t[Nel];
	89	if(rho>0.) fDensity = rho;
	90	else { // try to compute density form mixture.
	91	rho = 0.0;
	92	for(i=0;i<N;i++) if(nw[i]>0.0) {
	93	mix = (TGeoMixture*) (m->At(i));
	94	rho += nw[i]*(mix->GetDensity());
	95	} // end for i
	96	fDensity = rho;
	97	} // end if
	98	if(radlen>0.) fRadLen = radlen;
	99	else { // try to compute radiation form mixture.
	100	// From "Review of Particle Physics" Particle Data Group Section
	101	// 26.4.1 equation 26.21 (2002).
	102	radlen = 0.0;
	103	for(i=0;i<N;i++) if(nw[i]>0.0) {
	104	mix = (TGeoMixture*) (m->At(i));
	105	if(mix->GetRadLen()>0.0) rho += 1.0/(nw[i]*(mix->GetRadLen()));
	106	} // end for i
	107	fRadLen = 1.0/radlen;
	108	} // end if
	109	if(intleng>0.) fIntLen = intleng;
	110	else { // try to compute interaction form mixture.
	111	intleng = 0.0;
	112	for(i=0;i<N;i++) if(nw[i]>0.0) {
	113	mix = (TGeoMixture*) (m->At(i));
	114	if(mix->GetIntLen()>0.0) intleng += 1.0/(nw[i]*(mix->GetIntLen()));
	115	} // end for i
	116	fIntLen = 1.0/intleng;
	117	} // end if
	118	j = 0;
	119	for(z=1;z<110;z++){
	120	wel[z] /= tw;
	121	el[z] /= tw;
	122	Ael[z] /= tw;
	123	if(wel[z]>0.0) this->DefineElement(j++,Ael[z],el[z],wel[z]);
	124	} // end for i
	125	delete[] nw;
	126	}
	127	//======================================================================
	128	ClassImp(AliITSGeoCable)
	129	;
	130	AliITSGeoCable::AliITSGeoCable(){
	131	//
	132
	133	fRmin = fRmax = 0.0;
	134	fNs.SetXYZ(0.0,0.0,0.0);
	135	fNe.SetXYZ(0.0,0.0,0.0);
	136	fTubes = 0;
	137	fTranRot = 0;
	138	}
	139	//----------------------------------------------------------------------
	140	AliITSGeoCable::AliITSGeoCable(const char name,const TObjArray vect,
	141	const Double_t Rmin,const Double_t Rmax,
	142	const TVector3 ns,const TVector3 ne){
	143	//
	144	// Inputs:
	145	// char *name Name of this compound object
	146	// TObjArray *vect Array of TVector3's of points representing the
	147	// path of the cable
	148	// TVector3 ns=0 Normal vector representing the angle of the
	149	// starting surface, default perpendicular
	150	// TVector3 ne=0 Normal vector representing the angle of the
	151	// ending surface, default perpendicular
	152	// Outputs:
	153	// none.
	154	// Return:
	155	// A fully initilized and created AliITSGeoCable class.
	156	Char_t nam[500];
	157	Int_t i,n;
	158	Double_t s,th,ph;
	159	TVector3 x0,x1,x2,d,t,n0,n1;
	160
	161	fRmin = Rmin;
	162	fRmax = Rmax;
	163	fRmin = fRmax = 0.0;
	164	fNs.SetXYZ(0.0,0.0,-1.0);
	165	fNe.SetXYZ(0.0,0.0,1.0);
	166	n = vect->GetEntries();
	167	fTubes = new TObjArray(n-1);
	168	fTranRot = new TObjArray(n-1);
	169	fNs = ns*(1./ns.Mag());
	170	fNe = ne*(1./ne.Mag());
	171	//
	172	x0 = ((TVector3 )(vect->At(0)));
	173	n0 = ns;
	174	for(i=1;i<n;i++){
	175	x1 = ((TVector3 )(vect->At(i)));
	176	d = x1 - x0;
	177	if(i<n-1) {
	178	x2 = ((TVector3 )(vect->At(i+1)));
	179	n1 = d + (x2-x1);
	180	n1 *= 1./n1.Mag();
	181	}else{
	182	n1 = fNe;
	183	} // end if
	184	t = 0.5*(x1 + x0);
	185	th = TMath::ATan2(TMath::Sqrt(2.d.Mag2()-d.z()-2.d.Mag()*d.y()-
	186	2.d.Mag()d.x()),
	187	TMath::Sqrt(d.z()d.z()-2.d.Mag()*d.z()+d.Mag2()));
	188	th *= TMath::RadToDeg();
	189	ph = TMath::ATan2(d.y()-d.Mag(),d.x()-d.Mag());
	190	ph *= TMath::RadToDeg();
	191	sprintf(nam,"%sCombiTrans%dCable",name,i-1);
	192	fTranRot->AddAt(new TGeoCombiTrans(nam,t.x(),t.y(),t.z(),
	193	new TGeoRotation("",ph,th,0.0)),i-1);
	194	s = d.Mag();
	195	sprintf(nam,"%sPart%dCable",name,i-1);
	196	fTubes->AddAt( new TGeoCtub(nam,fRmin,fRmax,0.5*s,0.0,360.0,n0.x(),
	197	n0.y(),n0.z(),n1.x(),n1.y(),n1.z()),i-1);
	198	n0 = -n1;
	199	x0 = x1;
	200	} // end for i
	201	}
	202	//----------------------------------------------------------------------
	203	AliITSGeoCable::~AliITSGeoCable(){
	204	//
	205	// Inputs:
	206	// none.
	207	// Outputs:
	208	// none.
	209	// Return:
	210	// none.
	211	Int_t i;
	212
	213	if(fTubes){
	214	for(i=0;i<fTubes->GetEntries();i++)
	215	delete (TGeoCtub*)(fTubes->At(i));
	216	delete fTubes;
	217	} // end if
	218	fTubes = 0;
	219	if(fTranRot){
	220	for(i=0;i<fTranRot->GetEntries();i++)
	221	delete (TGeoCombiTrans*)(fTranRot->At(i));
	222	delete fTranRot;
	223	} // end if
	224	fTranRot = 0;
	225	}
	226	//----------------------------------------------------------------------