ComputeObsTime.C

Go to the documentation of this file.

/*
# Copyright (C) 2019 Gabriele Vedovato
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <https://www.gnu.org/licenses/>.
*/
 
//
// this macro is used to compute the coincidence and exclusive times starting from cat0 DQ
// it works for nifo<=3
// Consider for example the network with 3 detectoes H1, L1, V1, the following files are produced:
// 
// single detector cat0 observation time
// 
// -> H1_cat0.txt , L1_cat0.txt , V1_cat0.txt
// 
// double detectors cat0 coincident observation time
// 
// -> H1_L1_cat0.txt , H1_V1_cat0.txt , L1_V1_cat0.txt
// 
// triple detectors cat0 coincident observation time
// 
// -> H1_L1_V1_cat0.txt
// 
// double detectors cat0 exclusive observation time
// 
// -> H1_L1_xcat0.txt , H1_V1_xcat0.txt , L1_V1_xcat0.txt
// 
// single detector cat0 exclusive observation time
// 
// -> H1_xcat0.txt , L1_xcat0.txt , V1_xcat0.txt
//
// All informations and the total times for each network are dumped into a README.md file
//
 
#define DAY (24.*3600.)
 
#define MAX_IFOS  3
 
void ComputeObsTime(TString DQ_dir, TString OBSTIME_dir="") {
 
  // read cat0 file list
  vector<TString> ifname = CWB::Toolbox::getFileListFromDir(DQ_dir, ".txt","","_cat0",true);
  if(ifname.size()==0) {
    cout << "No cat0 files are present on the directory : " << DQ_dir << endl;
    exit(1);
  }
 
  // extract ifo list
  vector<TString> ifo;
  for(int j=0;j<ifname.size();j++) {
    //cout << j << " " << ifname[j] << endl;
    TString xifo = TString(gSystem->BaseName(ifname[j]))(0,2);
    ifo.push_back(xifo);
  }
 
  int nifo = ifo.size();
 
  if(nifo>MAX_IFOS) {
    cout << "Error, max ifos number is : " << MAX_IFOS << endl;
    exit(1);
  }
 
  // open README.md file
  ofstream out;
  if(CWB::Toolbox::isFileExisting("README.md")) {
    out.open("README.md",ios::app);   // append
  } else {
    out.open("README.md",ios::out);   // create
  }
  if(!out.good()) {cout << "ComputeObsTime.C - Error : Opening File : "
                        << "README.md" << endl;gSystem->Exit(1);}
 
  // compute single detector cat0 observation time
 
  out << "-------------------------------------------------------" << endl;
  out << "single detector cat0 observation time" << endl;
  out << "-------------------------------------------------------" << endl;
  out << endl; 
 
  vector<vector<waveSegment>> ifo1_seg(nifo);
  for(int i=0;i<nifo;i++) {
 
    ifo1_seg[i] = CWB::Toolbox::readSegments(ifname[i]);
 
    double ifo1_time=CWB::Toolbox::getTimeSegList(ifo1_seg[i]);
    out << ifo[i] << "_time = " << (int)ifo1_time << " sec - " << ifo1_time/DAY << " days" << endl;
 
    // dump single detector cat0 segments
    if(OBSTIME_dir!="") {
      TString ofname = TString::Format("%s/%s_cat0.txt",OBSTIME_dir.Data(),ifo[i].Data());
      out << "-> " << gSystem->BaseName(ofname) << endl << endl;
      CWB::Toolbox::dumpSegList(ifo1_seg[i],ofname); 
    } 
  }
  out << endl;
 
  // compute double detectors cat0 coincident observation time
 
  out << "-------------------------------------------------------" << endl;
  out << "double detectors cat0 coincident observation time" << endl;
  out << "-------------------------------------------------------" << endl;
  out << endl;
 
  vector<vector<vector<waveSegment>>> ifo2_seg(nifo);
  for(int i=0;i<nifo;i++) {
 
    ifo2_seg[i].resize(nifo);
 
    for(int j=i+1;j<nifo;j++) {
 
      ifo2_seg[i][j] = CWB::Toolbox::mergeSegLists(ifo1_seg[i],ifo1_seg[j]);
      double ifo2_time = CWB::Toolbox::getTimeSegList(ifo2_seg[i][j]);
      out << ifo[i] << "_" << ifo[j] << "_time = " << (int)ifo2_time << " sec - " << ifo2_time/DAY << " days" << endl;
 
      // dump double detector cat0 coincident segments
      if(OBSTIME_dir!="") {
        TString ofname = TString::Format("%s/%s_%s_cat0.txt",OBSTIME_dir.Data(),ifo[i].Data(),ifo[j].Data());
        out << "-> " << gSystem->BaseName(ofname) << endl << endl;
        CWB::Toolbox::dumpSegList(ifo2_seg[i][j],ofname); 
      } 
    }
  }
  out << endl;
 
  // compute triple detectors cat0 coincident observation time
 
  out << "-------------------------------------------------------" << endl;
  out << "triple detectors cat0 coincident observation time" << endl;
  out << "-------------------------------------------------------" << endl;
  out << endl;
 
  TString ifo3_name=ifo[0];
  vector<waveSegment> ifo3_seg = ifo1_seg[0];
  for(int i=1;i<nifo;i++) {
    ifo3_seg = CWB::Toolbox::mergeSegLists(ifo3_seg,ifo1_seg[i]);
    ifo3_name += "_"+ifo[i];
  }
  double ifo3_time = CWB::Toolbox::getTimeSegList(ifo3_seg);
  out << ifo3_name+"_time = " << (int)ifo3_time << " sec - " << ifo3_time/DAY << " days" << endl;
 
  if(OBSTIME_dir!="") {
    cout << endl;
    TString ofname = TString::Format("%s/%s_cat0.txt",OBSTIME_dir.Data(),ifo3_name.Data());
    out << "-> " << gSystem->BaseName(ofname) << endl << endl;
    CWB::Toolbox::dumpSegList(ifo3_seg,ofname); 
  }
  out << endl;
 
  vector<waveSegment> ifo3_iseg = CWB::Toolbox::invertSegments(ifo3_seg);  // compute inverse ifo3 segments
 
  // compute double detectors cat0 exclusive observation time
 
  out << "-------------------------------------------------------" << endl;
  out << "double detectors cat0 exclusive observation time" << endl;
  out << "-------------------------------------------------------" << endl;
  out << endl;
 
  vector<vector<vector<waveSegment>>> ifo2_xseg(nifo);
  vector<vector<vector<waveSegment>>> ifo2_ixseg(nifo);
  for(int i=0;i<nifo;i++) {
 
    ifo2_xseg[i].resize(nifo);
    ifo2_ixseg[i].resize(nifo);
 
    for(int j=i+1;j<nifo;j++) {
 
      ifo2_xseg[i][j] = CWB::Toolbox::mergeSegLists(ifo2_seg[i][j],ifo3_iseg);
      double ifo2_xtime = CWB::Toolbox::getTimeSegList(ifo2_xseg[i][j]);
      out << ifo[i] << "_" << ifo[j] << "_xtime = " << (int)ifo2_xtime << " sec - " << ifo2_xtime/DAY << " days" << endl;
 
      // dump double detector cat0 coincident segments
      if(OBSTIME_dir!="") {
        TString ofname = TString::Format("%s/%s_%s_xcat0.txt",OBSTIME_dir.Data(),ifo[i].Data(),ifo[j].Data());
        out << "-> " << gSystem->BaseName(ofname) << endl << endl;
        CWB::Toolbox::dumpSegList(ifo2_xseg[i][j],ofname); 
      } 
 
      ifo2_ixseg[i][j] = CWB::Toolbox::invertSegments(ifo2_xseg[i][j]);
    }
  }
  out << endl;
 
  // compute single detector cat0 exclusive observation time
 
  out << "-------------------------------------------------------" << endl;
  out << "single detector cat0 exclusive observation time" << endl;
  out << "-------------------------------------------------------" << endl;
  out << endl;
 
  for(int i=0;i<nifo;i++) for(int j=i+1;j<nifo;j++) ifo2_ixseg[j][i] = ifo2_ixseg[i][j];
 
  vector<vector<waveSegment>> ifo1_xseg(nifo);
  for(int i=0;i<nifo;i++) {
 
    ifo1_xseg[i] = CWB::Toolbox::mergeSegLists(ifo1_seg[i],ifo3_iseg);
 
    for(int j=0;j<nifo;j++) {
      if(j!=i) ifo1_xseg[i] = CWB::Toolbox::mergeSegLists(ifo1_xseg[i],ifo2_ixseg[i][j]);
    }
 
    double ifo1_xtime=CWB::Toolbox::getTimeSegList(ifo1_xseg[i]);
    out << ifo[i] << "_xtime = " << (int)ifo1_xtime << " sec - " << ifo1_xtime/DAY << " days" << endl;
 
    // dump single detector cat0 segments
    if(OBSTIME_dir!="") {
      TString ofname = TString::Format("%s/%s_xcat0.txt",OBSTIME_dir.Data(),ifo[i].Data());
      out << "-> " << gSystem->BaseName(ofname) << endl << endl;
      CWB::Toolbox::dumpSegList(ifo1_xseg[i],ofname); 
    } 
  }
  out << endl;
 
  out.close();
 
  exit(0);
}