Skip to content

Commit

Permalink
Revert "Merge pull request #424 from JeffersonLab/AddInsertToFCAL"
Browse files Browse the repository at this point in the history 
This reverts commit d8748b1, reversing
changes made to 79d8d16.
  • Loading branch information
@markito3
markito3 committed Aug 1, 2020
1 parent d8748b1 commit 6368820
Show file tree
Hide file tree
Showing 16 changed files with 279 additions and 447 deletions.
17 changes: 5 additions & 12 deletions src/libraries/FCAL/DFCALCluster.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -105,7 +105,7 @@ void DFCALCluster::resetClusterHits()
}

bool DFCALCluster::update( const userhits_t* const hitList,
double fcalFaceZ, const DFCALGeometry *fcalgeom )
double fcalFaceZ )
{

double energy = 0;
Expand Down Expand Up @@ -254,8 +254,7 @@ bool DFCALCluster::update( const userhits_t* const hitList,

for (int ih = 0; ih < hitList->nhits; ih++) {
shower_profile( hitList, ih,fEallowed[ih],fEexpected[ih],
fcalFaceZ+0.5*DFCALGeometry::blockLength(),
fcalgeom);
fcalFaceZ+0.5*DFCALGeometry::blockLength());
}
}

Expand All @@ -265,7 +264,7 @@ bool DFCALCluster::update( const userhits_t* const hitList,
void DFCALCluster::shower_profile( const userhits_t* const hitList,
const int ihit,
double& Eallowed, double& Eexpected,
double fcalMidplaneZ, const DFCALGeometry *fcalgeom) const
double fcalMidplaneZ) const
{

//std::cout << " Run profile for hit " << ihit;
Expand All @@ -274,12 +273,6 @@ void DFCALCluster::shower_profile( const userhits_t* const hitList,
return;
double x = hitList->hit[ihit].x;
double y = hitList->hit[ihit].y;
double moliere_radius=MOLIERE_RADIUS;
double min_dist=fcalgeom->blockSize();
if (fabs(x)<fcalgeom->insertSize() && fabs(y)<fcalgeom->insertSize()){
moliere_radius=PbWO4_MOLIERE_RADIUS;
min_dist=fcalgeom->insertBlockSize();
}
double dist = sqrt(SQR(x - fCentroid.x()) + SQR(y - fCentroid.y()));
if (dist > MAX_SHOWER_RADIUS)
return;
Expand All @@ -289,7 +282,7 @@ void DFCALCluster::shower_profile( const userhits_t* const hitList,
double v0 = 0;
double u = x*cos(phi) + y*sin(phi);
double v =-x*sin(phi) + y*cos(phi);
double vVar = SQR(moliere_radius);
double vVar = SQR(MOLIERE_RADIUS);
double uVar = vVar+SQR(SQR(8*theta));
double vTail = 4.5+0.9*log(fEnergy+0.05);
double uTail = vTail+SQR(10*theta);
Expand All @@ -298,7 +291,7 @@ void DFCALCluster::shower_profile( const userhits_t* const hitList,
Eexpected = fEnergy*core;
Eallowed = 2*fEmax*core + (0.2+0.5*log(fEmax+1.))*tail;

if ((dist <= min_dist) && (Eallowed < fEmax) ) {
if ((dist <= 4.) && (Eallowed < fEmax) ) {
std::cerr << "Warning: FCAL cluster Eallowed value out of range!\n";
Eallowed = fEmax;
}
Expand Down
9 changes: 2 additions & 7 deletions src/libraries/FCAL/DFCALCluster.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -16,11 +16,8 @@ using namespace std;
#include <JANA/JFactory.h>
using namespace jana;

#include <FCAL/DFCALGeometry.h>

#define FCAL_USER_HITS_MAX 2800
#define MOLIERE_RADIUS 3.696
#define PbWO4_MOLIERE_RADIUS 2.2
#define MAX_SHOWER_RADIUS 25

class DFCALCluster : public JObject {
Expand Down Expand Up @@ -81,8 +78,7 @@ class DFCALCluster : public JObject {
int getHits() const; // get number of hits owned by a cluster
int addHit(const int ihit, const double frac);
void resetClusterHits();
bool update( const userhits_t* const hitList, double fcalFaceZ,
const DFCALGeometry *fcalgeom );
bool update( const userhits_t* const hitList, double fcalFaceZ );

// get hits that form a cluster after clustering is finished
inline const vector<DFCALClusterHit_t> GetHits() const { return my_hits; }
Expand All @@ -101,8 +97,7 @@ class DFCALCluster : public JObject {
void shower_profile( const userhits_t* const hitList,
const int ihit,
double& Eallowed, double& Eexpected,
double fcalMidplaneZ,
const DFCALGeometry *fcalgeom ) const ;
double fcalMidplaneZ ) const ;

// internal parsers of properties for a hit belonging to a cluster
oid_t getHitID( const userhits_t* const hitList, const int ihit) const;
Expand Down
11 changes: 6 additions & 5 deletions src/libraries/FCAL/DFCALCluster_factory.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -96,8 +96,9 @@ jerror_t DFCALCluster_factory::evnt(JEventLoop *eventLoop, uint64_t eventnumber)
// try clusterizing if more than 250 hits in FCAL
if(fcalhits.size() > MAX_HITS_FOR_CLUSTERING) return NOERROR;

const DFCALGeometry* fcalGeom=NULL;
eventLoop->GetSingle(fcalGeom);
vector<const DFCALGeometry*> geomVec;
eventLoop->Get(geomVec);
const DFCALGeometry& fcalGeom = *(geomVec[0]);

// Sort hits by energy
sort(fcalhits.begin(), fcalhits.end(), FCALHitsSort_C);
Expand All @@ -112,7 +113,7 @@ jerror_t DFCALCluster_factory::evnt(JEventLoop *eventLoop, uint64_t eventnumber)
hit != fcalhits.end(); hit++ ) {
if ( (**hit).E < 1e-6 ) continue;
hits->hit[nhits].id = (**hit).id;
hits->hit[nhits].ch = fcalGeom->channel( (**hit).row, (**hit).column );
hits->hit[nhits].ch = fcalGeom.channel( (**hit).row, (**hit).column );
hits->hit[nhits].x = (**hit).x;
hits->hit[nhits].y = (**hit).y;
hits->hit[nhits].E = (**hit).E;
Expand Down Expand Up @@ -146,7 +147,7 @@ jerror_t DFCALCluster_factory::evnt(JEventLoop *eventLoop, uint64_t eventnumber)
bool something_changed = false;
for ( unsigned int c = 0; c < clusterCount; c++ ) {
//cout << " --------- Update iteration " << iter << endl;
something_changed |= clusterList[c]->update( hits, fcalFaceZ_TargetIsZeq0, fcalGeom );
something_changed |= clusterList[c]->update( hits, fcalFaceZ_TargetIsZeq0 );
}
if (something_changed) {
for ( unsigned int c = 0; c < clusterCount; c++ ) {
Expand Down Expand Up @@ -177,7 +178,7 @@ jerror_t DFCALCluster_factory::evnt(JEventLoop *eventLoop, uint64_t eventnumber)
clusterList[clusterCount] = new DFCALCluster( hits->nhits );
hitUsed[ih] = -1;
clusterList[clusterCount]->addHit(ih,1.);
clusterList[clusterCount]->update( hits, fcalFaceZ_TargetIsZeq0, fcalGeom );
clusterList[clusterCount]->update( hits, fcalFaceZ_TargetIsZeq0 );
++clusterCount;
}
else if (iter > 0) {
Expand Down
170 changes: 49 additions & 121 deletions src/libraries/FCAL/DFCALGeometry.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -15,127 +15,70 @@ using namespace std;
//---------------------------------
// DFCALGeometry (Constructor)
//---------------------------------
DFCALGeometry::DFCALGeometry(const DGeometry *geom){
double innerRadius = ( kBeamHoleSize - 1 ) / 2. * blockSize() * sqrt(2.);

// inflate the innner radius by 1% to for "safe" comparison
innerRadius *= 1.01;

// Check for presence of PbWO4 insert
int insert_row_size=0;
geom->Get("//composition[@name='LeadTungstateFullRow']/mposX[@volume='LTBLwrapped']/@ncopy",insert_row_size);
m_insertSize=insertBlockSize()*double(insert_row_size/2);

geom->GetFCALPosition(m_FCALdX,m_FCALdY,m_FCALfront);
DVector2 XY0(m_FCALdX,m_FCALdY);

vector<double>block;
geom->Get("//box[@name='LGBL']/@X_Y_Z",block);
double back=m_FCALfront+block[2];
geom->Get("//box[@name='LTB1']/@X_Y_Z",block);
m_insertFront=0.5*(back+m_FCALfront-block[2]);


// Initilize the list of active blocks to false, to be adjusted for the
// actual geometry below.
for( int row = 0; row < 2*kBlocksTall; row++ ){
for( int col = 0; col < 2*kBlocksWide; col++ ){
m_activeBlock[row][col] = false;
}
}
DFCALGeometry::DFCALGeometry() :
m_numActiveBlocks( 0 )
{
double innerRadius = ( kBeamHoleSize - 1 ) / 2. * blockSize() * sqrt(2.);

// Now fill in the data for the actual geometry
unsigned int ch=0;
for( int row = 0; row < kBlocksTall; row++ ){
for( int col = 0; col < kBlocksWide; col++ ){
// inflate the innner radius by 1% to for "safe" comparison
innerRadius *= 1.01;

for( int row = 0; row < kBlocksTall; row++ ){
for( int col = 0; col < kBlocksWide; col++ ){

// transform to beam axis
m_positionOnFace[row][col] =
DVector2( ( col - kMidBlock ) * blockSize(),
( row - kMidBlock ) * blockSize());

double thisRadius = m_positionOnFace[row][col].Mod();
// transform to beam axis
m_positionOnFace[row][col] =
DVector2( ( col - kMidBlock ) * blockSize(),
( row - kMidBlock ) * blockSize() );

if( ( thisRadius < radius() ) && ( thisRadius > innerRadius )
){
// Carve out region for insert, if present. For back compatibility
// we set these blocks as inactive but maintain the length of the array
if (fabs(m_positionOnFace[row][col].X())>m_insertSize ||
fabs(m_positionOnFace[row][col].Y())>m_insertSize){
m_activeBlock[row][col] = true;
m_positionOnFace[row][col]+=XY0; // add FCAL offsets
}
// build the "channel map"
m_channelNumber[row][col] = ch;
m_row[ch] = row;
m_column[ch] = col;

ch++;
}
}
}

if (insert_row_size>0){
m_insertMidBlock=(insert_row_size-1)/2;
for( int row = 0; row < insert_row_size; row++ ){
for( int col = 0; col < insert_row_size; col++ ){

// transform to beam axis
int row_index=row+kBlocksTall;
int col_index=col+kBlocksWide;
m_positionOnFace[row_index][col_index] =
DVector2( ( col - m_insertMidBlock -0.5) * insertBlockSize(),
( row - m_insertMidBlock -0.5) * insertBlockSize() );

if( fabs(m_positionOnFace[row_index][col_index].X())>insertBlockSize()
|| fabs(m_positionOnFace[row_index][col_index].Y())>insertBlockSize()
){
m_activeBlock[row_index][col_index] = true;
m_positionOnFace[row_index][col_index]+=XY0; // add FCAL offsets

// build the "channel map"
m_channelNumber[row_index][col_index] = ch;
m_row[ch] = row_index;
m_column[ch] = col_index;

ch++;
double thisRadius = m_positionOnFace[row][col].Mod();

if( ( thisRadius < radius() ) && ( thisRadius > innerRadius ) ){

m_activeBlock[row][col] = true;

// build the "channel map"
m_channelNumber[row][col] = m_numActiveBlocks;
m_row[m_numActiveBlocks] = row;
m_column[m_numActiveBlocks] = col;

m_numActiveBlocks++;
}
else{

m_activeBlock[row][col] = false;
}
}
}
}
}
}
m_numChannels=ch;
}

bool
DFCALGeometry::isBlockActive( int row, int column ) const
{
if (row>=100&&column>=100){
row-=100-kBlocksTall;
column-=100-kBlocksWide;
}
// I'm inserting these lines to effectively disable the
// two assert calls below. They are causing all programs
// (hd_dump, hdview) to exit, even when I'm not interested
// in the FCAL. This does not fix the underlying problem
// of why we're getting invalid row/column values.
// 12/13/05 DL
if( row < 0 || row >= kBlocksTall )return false;
if( column < 0 || column >= kBlocksWide )return false;

// assert( row >= 0 && row < kBlocksTall );
// assert( column >= 0 && column < kBlocksWide );

return m_activeBlock[row][column];
}

int
DFCALGeometry::row( float y, bool in_insert ) const
DFCALGeometry::row( float y ) const
{
y-=m_FCALdY;

if (in_insert){
return kBlocksTall+static_cast<int>( y / insertBlockSize() + m_insertMidBlock + 0.5);
}
return static_cast<int>( y / blockSize() + kMidBlock + 0.5);
}

int
DFCALGeometry::column( float x, bool in_insert ) const
DFCALGeometry::column( float x ) const
{
x-=m_FCALdX;

if (in_insert){
return kBlocksWide+static_cast<int>( x / insertBlockSize() + m_insertMidBlock + 0.5);
}
return static_cast<int>( x / blockSize() + kMidBlock + 0.5);
}

Expand All @@ -144,28 +87,21 @@ DFCALGeometry::positionOnFace( int row, int column ) const
{
// assert( row >= 0 && row < kBlocksTall );
// assert( column >= 0 && column < kBlocksWide );
// Check for insert blocks
if (row>=100&&column>=100){
row-=100-kBlocksTall;
column-=100-kBlocksWide;
}
return m_positionOnFace[row][column];

return m_positionOnFace[row][column];
}

DVector2
DFCALGeometry::positionOnFace( int channel ) const
{
return positionOnFace( m_row[channel], m_column[channel] );
assert( channel >= 0 && channel < m_numActiveBlocks );

return positionOnFace( m_row[channel], m_column[channel] );
}

int
DFCALGeometry::channel( int row, int column ) const
{
// Check for insert blocks
if (row>=100&&column>=100){
row-=100-kBlocksTall;
column-=100-kBlocksWide;
}
if( isBlockActive( row, column ) ){

return m_channelNumber[row][column];
Expand All @@ -177,11 +113,3 @@ DFCALGeometry::channel( int row, int column ) const
return -1;
}
}

bool DFCALGeometry::inInsert(int channel) const{
if (fabs(positionOnFace(channel).X()-m_FCALdX)<m_insertSize
&& fabs(positionOnFace(channel).Y()-m_FCALdY)<m_insertSize){
return true;
}
return false;
}
Loading

0 comments on commit 6368820

Please sign in to comment.