quicklook/QLflightS4_ND/S4_QL.cpp


/*
* S4 Quick Look
* Author Marcelli-Malvezzi
* Version 1.00 - March 2006
* 
* Description - The aim of S4 QL software is to monitor the bahaviour of this detector.
*               It creates three canvases: the first one contains the histogram of S4 data, the second one 
*               is relative to the time behaviour of collected data while the last shows S4 rate 
*               (information from trigger Packet).
*               See documentation for a more detailed description of the output.
*
*
* Parameters:
*               TSTring base    - the path to the root directory for the specific Pamela unpack session
*                                 There is no default value, without this input the program will not run
*               TString outDir  - the path where to save the output image (Default = ./) 
*               TString format  - the format which will be used fo rsave the produced images (Default = "jpg") 
*               Float_t DeltaT  - the time interval in minute for calculation of average S4 data for minute, 
*                                 see S4_QL_2 plot (Default = 1 minute)
*
*
* Version 1.1  - June 2006
* Fixed bugs: the vector "trcss" was inizilized to have dimension 10, but for files with large number of events this
*             is not true; now this vector is inizialized at 100.
*
*             the threshold at which S4 is set and the trigger configuration can change in the file; all these changes are reported
*             in a pad
*
*             for a large namber of events is not possible to have vectors, so all graphs have been converted in histograms
*
* Version 2.0 - September 2006
* Fixed bugs: changed the number of bin in the "rate" histograms and the size of the string vector trcss
*******/


#include <iostream>
#include <fstream>
#include <sstream>
#include <math.h>
#include "TStyle.h"
#include "TFile.h"
#include "TList.h"
#include "TTree.h"
#include "TLatex.h"
#include "TObjString.h"
#include "TCanvas.h"
#include "TGraph.h"
#include "TH1F.h"
#include "TF1.h"
#include "TH2F.h"
#include "TF2.h"
#include "TGaxis.h"
#include "TString.h"
#include "TPaveText.h"
#include "EventHeader.h"
#include "PscuHeader.h"
#include "TMultiGraph.h"
#include "physics/S4/S4Event.h"
#include "varDump/VarDumpEvent.h"
#include "varDump/VarDumpRecord.h"
#include "physics/trigger/TriggerEvent.h"

using namespace std;

void S4_QL(TString base, TString outDir, TString format, ULong_t DeltaT){  //DeltaT in minute
  
  //------ Variables initialization ---------/
  Int_t tmpSize;
  ULong_t mintime, maxtime;
  Int_t   adcmax;
  Int_t   j=0;
  Int_t   S4_TRHadc;
  Int_t   S4_TRHmip;
  TString filename = ((TObjString*)base.Tokenize('/')->Last())->GetString();
  filename = ((TObjString*)filename.Tokenize('.')->First())->GetString();
  char *trc;
  TString trcs,ciao;
  TString str;
  TString trcsstot[31];
  TString trcsstot2;
  
  TString trgconf[31]={"TOF1","TOF2","TOF3","TOF4","TOF5","TOF6","TOF7","S4","CALO","CALO-S4","TOF1-S4","TOF2-S4","TOF3-S4","TOF4-S4","TOF5-S4","TOF6-S4","TOF7-S4","TOF1-CALO","TOF2-CALO","TOF3-CALO","TOF4-CALO","TOF5-CALO","TOF6-CALO","TOF7-CALO","TOF1-CALO-S4","TOF2-CALO-S4","TOF3-CALO-S4","TOF4-CALO-S4","TOF5-CALO-S4","TOF6-CALO-S4","TOF7-CALO-S4",};
  stringstream oss, oss1, oss2, oss3, s4soglia, buffer, conftrig, noentries;
  ULong_t lastime, firstime,obmin,obmax,obt;
  Int_t vardumpentries = 0;
  Int_t channels = 4096;
  char S4_TRH[10];
  char S4_TRH2[10];
  Int_t trigconf = 0;
  Long64_t nevents;
  string title;
    
  double s4rate;
  //------to open headerfile, s4file, Trigfile and vardumpfile---------------------------/
  TFile  *file  =new TFile(base.Data()) ;
  if (!file){
    printf("No such file in the directory has been found");
    return;
  }
  
  TTree *VarDumpTr = (TTree*)file->Get("VarDump");
  TTree *PhysicsTr = (TTree*)file->Get("Physics");
  
  TBranch *S4Br = PhysicsTr->GetBranch("S4");
  TBranch *TriggerBr = PhysicsTr->GetBranch("Trigger");
  TBranch *headBr = PhysicsTr->GetBranch("Header");
  TBranch *VarDumpBr = VarDumpTr->GetBranch("VarDump");
  TBranch *headVarDumpBr = VarDumpTr->GetBranch("Header");
  
  pamela::S4::S4Event  *s4e  = 0;
  pamela::EventHeader  *eh   = 0;
  pamela::PscuHeader   *ph   = 0;
  pamela::trigger::TriggerEvent *trige = 0;
  pamela::VarDumpEvent  *vde = 0; 
  pamela::VarDumpRecord *vdr = 0;
 
  PhysicsTr->SetBranchAddress("S4", &s4e);
  PhysicsTr->SetBranchAddress("Header", &eh);
  PhysicsTr->SetBranchAddress("Trigger", &trige);
  VarDumpTr->SetBranchAddress("VarDump", &vde);

  nevents = S4Br->GetEntries();

  //----------- If nevents < = 0 ---------------------------------/
  if (nevents<=0) {
    printf("nevents = %i \n", nevents);
    printf(" \n");
    
    TCanvas *canvas4 = new TCanvas("No entries", "No entries ", 400, 200);
    canvas4->SetFillColor(10);
    canvas4->cd();

    TLatex *l = new TLatex();
    l->SetTextAlign(12);
    l->SetTextSize(0.15);
    l->SetTextColor(2);
    noentries.str("");
    noentries<< "S4_QL:";
    l->DrawLatex(0.05, 0.7, noentries.str().c_str()); 
    noentries.str("");
    noentries<< "No entries for this files";
    l->DrawLatex(0.05, 0.5, noentries.str().c_str()); 
    
    if (outDir == "./") {
      oss.str("");
      oss << filename.Data() << "_S4_QL." << format.Data();
    } else {
      oss.str("");
      oss << outDir.Data() << filename.Data() << "_S4_QL." << format.Data();
    }
    
    canvas4->Update();
    canvas4->SaveAs(oss.str().c_str()); 
    
    return;
  }
  const Int_t ev = nevents;
  TString trcss[ev]; 
  //-------to set trigger configuration---------------------------------//
  for (Int_t ll = 0; ll < nevents; ll++){
    TriggerBr->GetEntry(ll);
    trigconf = trige->trigconf;
    trc = 0;
    if ( trigconf & (1<<0) ) trc = "TOF1"; 
    if ( trigconf & (1<<1) ) { 
      if (trc==0) trc= "TOF2";
      else trc = Form("%s-TOF2",trc);
    }
    if ( trigconf & (1<<2) ) { 
      if (trc==0) trc= "TOF3";
      else trc = Form("%s-TOF3",trc);
    } 
    if ( trigconf & (1<<3) ) { 
      if (trc==0) trc= "TOF4";
      else trc = Form("%s-TOF4",trc);
    }  
    if ( trigconf & (1<<4) ) { 
      if (trc==0) trc= "TOF5";
      else trc = Form("%s-TOF5",trc);
    }  
    if ( trigconf & (1<<5) ) { 
      if (trc==0) trc= "TOF6";
      else trc = Form("%s-TOF6",trc);
    }
    if ( trigconf & (1<<6) ) { 
      if (trc==0) trc= "TOF7";
      else trc = Form("%s-TOF7",trc);
    } 
    if ( trigconf & (1<<7) ) { 
      if (trc==0) trc= "S4";
      else trc = Form("%s-S4",trc);
    } 
   if ( trigconf & (1<<8) ) { 
     if (trc==0) trc= "CALO";
     else trc = Form("%s-CALO",trc);
   } 
   if ( trigconf & (1<<9) ) { 
     if (trc==0) trc= "CALIB_ON";
     else trc = Form("%s-CALIB_ON",trc);
   }  
   trcs = "";
   trcs = trc;
   if (ll==0){ 
     trcss[0]=trcs; 
     j=j+1;
   } 
   if (trcs!=trcss[j-1] && ll>0){ 
     trcss[j]=trcs;
     j=j+1; 
   } 
   //cout<<i<<"\n";
   //cout<<nevents<<"\n\n";
  }
  //--------------------------//
  for(Int_t k=0; k<31; k++){
    for(Int_t p=0; p<j ; p++){
      if(trgconf[k]==trcss[p]){
        trcsstot[k]=trgconf[k];
      }
    }
  }
  for(Int_t k=0; k<31; k++){
    if(trcsstot[k]!= "")
      trcsstot2=trcsstot2+"/"+ trcsstot[k];
  }
  
  //---------to search S4 threshold and convert it into char (hexadecimal value)--------//
  vardumpentries =  VarDumpBr->GetEntries();
  if(vardumpentries==0){
    cout<<"TTree VarDump: Entries = 0"<<"\n"<<"         S4 Threshold not defined"<<"\n";
    str= " Not defined";
  }
  else{
    for(Int_t k=0; k< vardumpentries; k++){
      VarDumpBr->GetEntry(k);
      vdr = (pamela::VarDumpRecord*)vde->Records->At(118);
      S4_TRHadc=((vdr->VAR_VALUE)/256);      //trh (number of ADC channels) 
      S4_TRHmip=S4_TRHadc*2;  

      //-----to convert decimal value to TString of S4 threshold value-------//
      if(k==0){
        sprintf(S4_TRH, "%d" , S4_TRHmip);
        str= S4_TRH;
      }else{
        sprintf(S4_TRH2, "%d" , S4_TRHmip);
        if(!strcmp(S4_TRH, S4_TRH2)) continue;
        str= str+"/"+S4_TRH2;   
      }
      sprintf(S4_TRH, "%d" , S4_TRHmip);
    }
    } 

  //********************** First Histogram ************************************//
  TH1F *h1 = new TH1F("All events", "S4 distribution for file: "+filename, channels, 1, channels);
  TH1F *h2 = new TH1F("S4 Triggered events", "S4 distribution for file: "+filename, channels, 1, channels);

  //********************** Second Histogram ************************************//
  headBr->GetEntry(0);
  ph = eh->GetPscuHeader();
  firstime = ph->GetOrbitalTime();

  headBr->GetEntry(nevents);
  ph = eh->GetPscuHeader();
  lastime = ph->GetOrbitalTime();
  headBr->GetEntry(nevents-100);
  ph = eh->GetPscuHeader();
  obt=ph->GetOrbitalTime();
  if(lastime-obt> 100000)nevents=nevents-100;
  int iii=0,MaxR=0;
  while(iii < nevents){
    headBr->GetEntry(iii);
    ph = eh->GetPscuHeader();
    obt=ph->GetOrbitalTime();
    if(obt<=firstime) firstime=ph->GetOrbitalTime();
    if(obt>=lastime) lastime=ph->GetOrbitalTime();
    TriggerBr->GetEntry(iii);
    if(trige->s4calcount[0]>MaxR)MaxR=trige->s4calcount[0];
    iii++;
  }
  MaxR=2*(1.05*MaxR/2);

  obmin=firstime;
  obmax=lastime;   
  const ULong_t nint=(((lastime-firstime)/(DeltaT*60000)));
  const Int_t size = (Int_t)(nint+1); 

  TH1F *Allev = new TH1F("Mean signal from S4-all events", filename+": Mean signal from S4 (all triggered events)", size, obmin, obmax);
  TH1F *Alltime = new TH1F("Mean signal from S4-all events", filename+": Mean signal from S4 (all triggered events)", size, obmin, obmax);
  TH1F *S4ev = new TH1F("Mean signal from S4-triggered events", filename+": Mean signal from S4 (only s4 triggered events)", size, obmin, obmax);
  TH1F *S4time = new TH1F("Mean signal from S4-triggered events", filename+": Mean signal from S4 (only s4 triggered events)", size, obmin, obmax);  

  //********************** Third Histogram ************************************//
  title="";
  title=filename+": S4 rate from Trigger Packet";

  const ULong_t nint2=(lastime-firstime)/10000;
  const Int_t size2 = (Int_t)((nint2+1)); // one bin every 10 second
  const ULong_t nint3=(lastime-firstime)/1000;
  const Int_t size3 = (Int_t)((nint3+1)); // one bin every 1 second
  const Int_t Rbin = MaxR/2;
  TH2I *rate= new TH2I(title.c_str(), title.c_str(), size2,obmin,obmax,Rbin,0,MaxR);
  TH2I *rateline= new TH2I(filename+". S4 rate from Trigger Packet: mean value over 100 events", filename+". S4 rate from Trigger Packet: mean value over 100 events",size3,obmin,obmax,Rbin,0,MaxR);

  //------------------------------------------------------------------------------------------------------------------
  //------- fill histograms ---------//
  Int_t n=0, p=0;
  for (Int_t i = 0; i < nevents; i++){
    TriggerBr->GetEntry(i);
    headBr->GetEntry(i);
    S4Br->GetEntry(i);
    ph = eh->GetPscuHeader();
    if (s4e->unpackError == 1 && (s4e->S4_DATA)==0)  continue;   
    obt=ph->GetOrbitalTime();
    rate->Fill(obt,trige->s4calcount[0]);
    h1->Fill(s4e->S4_DATA);  
    Allev->Fill(obt,s4e->S4_DATA);
    Alltime->Fill(obt);
    if ((trige->patterntrig[0] == 0)&&(trige->patterntrig[1] != 0)&&(trige->patterntrig[2] == 0)&&(trige->patterntrig[3] == 0)&&(trige->patterntrig[4] == 0)&&(trige->patterntrig[5] == 0)){
      h2->Fill(s4e->S4_DATA);
      S4ev->Fill(obt,s4e->S4_DATA);
      S4time->Fill(obt);
      p=p+1;
    }
  }
  Int_t kk=0;
  ULong_t Aobt=0;
  while (kk < nevents){
    obt=0;
    s4rate=0;
    Aobt=0;
    for(Int_t jj=kk; jj< (kk+100); jj++){
      TriggerBr->GetEntry(jj);
      headBr->GetEntry(jj);
      obt=ph->GetOrbitalTime();  
      s4rate= s4rate+(trige->s4calcount[0]);
      Aobt=Aobt-obmin+obt;
        //*(trige->s4calcount[0])/1000;
    }
    Aobt=Aobt/100;
    Aobt=Aobt+obmin;
    rateline->Fill(Aobt,s4rate/100);
    //rateline->Fill(ph->GetOrbitalTime(),(s4rate/100));
    kk=kk+100;
  }

  //****************************** Canvases *******************************//
  //------------------- First Canvas --------------------------------//
  TCanvas *canvas1 = new TCanvas("S4_QL_1", "S4 HISTO ", 1280, 1024);
  canvas1->SetFillColor(10);
  canvas1->Divide(1,2);

  TPad *all= new TPad ("","", 0, 0, 1, 1);
  all->SetFillColor(10);
  TPad *s4 = new TPad ("s4","s4", 0, 1, 1, 0); 
  s4->SetFillColor(10);
  TPad *pad = new TPad("pad","pad", .30,.80,.78,.90);
  pad->SetFillColor(10);
  TPad *pad2 = new TPad("pad2","pad2", .30,.80,.78,.90); 
  pad2->SetFillColor(10);
  TLatex *l = new TLatex();
  l->SetTextAlign(12);
  l->SetTextSize(0.35);
  l->SetTextColor(8);
  s4soglia.str("");
  s4soglia << "S4_THRESHOLD: "<<str.Data()<<" MIP";
  TLatex *ll = new TLatex();
  ll->SetTextAlign(12);
  ll->SetTextSize(0.35);
  conftrig.str("");
  conftrig<<"TRIGGER CONF: "<<trcsstot2.Data();
  if(adcmax<=100) adcmax=120;

  canvas1->cd(1);
  all->Draw();
  all->cd();
  h1->SetLineColor(kBlack);
  h1->SetFillColor(kRed);
  h1->GetXaxis()->SetTitle("ADC");
  h1->GetXaxis()->CenterTitle();
  h1->GetYaxis()->SetTitle("Number of Events");
  h1->GetYaxis()->CenterTitle(); 
  h1->Draw();
  pad->Draw();
  pad->cd();
  l->DrawLatex(0.05, 0.65, s4soglia.str().c_str());
  ll->DrawLatex(0.05, 0.15, conftrig.str().c_str());


  canvas1->cd(2);
  s4->Draw();
  s4->cd();
  h2->SetLineColor(kBlack);
  h2->SetFillColor(2);
  h2->GetXaxis()->SetTitle("ADC");
  h2->GetXaxis()->CenterTitle();
  h2->GetYaxis()->SetTitle("Number of Events");
  h2->GetYaxis()->CenterTitle(); 
  h2->Draw();
  pad2->Draw();
  pad2->cd();
  l->DrawLatex(0.05, 0.65, s4soglia.str().c_str()); 
  ll->DrawLatex(0.05, 0.15, conftrig.str().c_str());

  if (h1->GetMaximum()>0){
    all->SetLogy();
    all->SetLogx();
  } 
  if (h2->GetMaximum()>0){
    s4->SetLogy();
    s4->SetLogx();
  }
  canvas1->Update();
 

  if (outDir == "./") {
    oss1.str("");
    oss1 << filename.Data() << "_S4_QL_1." << format.Data();
  } else {
    oss1.str("");
    oss1 << outDir.Data() << filename.Data() << "_S4_QL_1." << format.Data(); 
  }
  canvas1->SaveAs(oss1.str().c_str()); 
 
  //------------------- Second Canvas --------------------------------//
  if (nint==0){
    cout<<"Number of Time Interval = 0"<<"\n"<<"      Set another value for input parameter DeltaTevtime "<<"\n"<<"RETURN"<<"\n";
    
    TCanvas *canvas2 = new TCanvas("Number of time interval=0", "Number of time interval=0", 900, 300);
    canvas2->SetFillColor(10);
    canvas2->cd();

    TLatex *l = new TLatex();
    l->SetTextAlign(12);
    l->SetTextSize(0.10);
    l->SetTextColor(2);
    noentries.str("");
    noentries<< "S4_QL: Time evolution";
    l->DrawLatex(0.05, 0.7, noentries.str().c_str()); 
    noentries.str("");
    noentries<< "Number of Time Interval = 0";
    l->SetTextColor(1);
    l->DrawLatex(0.05, 0.5, noentries.str().c_str()); 
    noentries.str("      Set another value for input parameter DeltaTevtime");
    noentries<< "";
    l->DrawLatex(0.05, 0.3, noentries.str().c_str()); 

    if (outDir == "./") {
      oss.str("");
      oss << filename.Data() << "_S4_QL_2." << format.Data();
    } else {
      oss.str("");
      oss << outDir.Data() << filename.Data() << "_S4_QL_2." << format.Data();
    }
    canvas2->Update();
    canvas2->SaveAs(oss.str().c_str()); 

    return;
  }
  TCanvas *canvas2 = new TCanvas("S4_QL_2", "S4 - Time evolution", 1280, 1024);
  canvas2->SetFillColor(10);
  canvas2->Divide(1,2);
  
  TPad *pad3 = new TPad("pad3","pad3", .7,.1,.88,.30);
  pad3->SetFillColor(10);
  TPad *pad4 = new TPad("pad4","pad4", .7,.1,.88,.30); 
  pad4->SetFillColor(10);

  TLatex *l3 = new TLatex();
  l3->SetTextAlign(12);
  l3->SetTextSize(.3);
  l3->SetTextColor(2);
  buffer.str("");
  buffer<<"DT :  "<<DeltaT<<" min";
  
  canvas2->cd(1);
  Allev->SetStats(kFALSE);
  Allev->SetMarkerColor(2);
  Allev->SetMarkerSize(.5);
  Allev->SetMarkerStyle(21);
  Allev->GetXaxis()->SetTitle("Time ( OBT in ms )");
  Allev->GetXaxis()->CenterTitle();
  Allev->GetYaxis()->SetTitle("Mean signal ( ADC )");
  Allev->GetYaxis()->CenterTitle();
  //Allev->SetMinimum(20);
  Allev->Divide(Alltime);
  Allev->Draw("p");
  pad3->Draw();
  pad3->cd();
  l3->DrawLatex(0.05, 0.65, buffer.str().c_str());

  
  canvas2->cd(2);
  if(p != 0){
  S4ev->SetStats(kFALSE);
  S4ev->SetMarkerColor(2);
  S4ev->SetMarkerSize(.5);
  S4ev->SetMarkerStyle(21);
  S4ev->GetXaxis()->SetTitle("Time ( OBT in ms )");
  S4ev->GetXaxis()->CenterTitle();
  S4ev->GetYaxis()->SetTitle("Mean signal ( ADC )");
  S4ev->GetYaxis()->CenterTitle();
  S4ev->SetMinimum(20);
  S4ev->Divide(S4time);
  S4ev->Draw("p");
  pad4->Draw();
  pad4->cd();
  l3->DrawLatex(0.05, 0.65, buffer.str().c_str());
  }else{
    TLatex *l = new TLatex();
    l->SetTextAlign(12);
    l->SetTextSize(0.07);
    l->SetTextColor(2);
    noentries.str("");
    noentries<< "S4_QL: Time evolution";
    l->DrawLatex(0.06, 0.6, noentries.str().c_str()); 
    noentries.str("");
    noentries<< "No S4 triggered events.";
    l->SetTextColor(1);
    l->DrawLatex(0.06, 0.5, noentries.str().c_str()); 
  }

  if (outDir == "./") {
    oss2.str("");
    oss2 << filename.Data() << "_S4_QL_2." << format.Data();
  } else {
    oss2.str("");
    oss2 << outDir.Data() << filename.Data() << "_S4_QL_2." << format.Data();
  }
  canvas2->SaveAs(oss2.str().c_str());

  //------------------- Third Canvas --------------------------------//
  TCanvas *canvas3 = new TCanvas("S4 - Rate", "S4 - Rate", 1280, 1024);
  canvas3->SetFillColor(10);
  canvas3->Divide(1,2);

  canvas3->cd(1);
  rate->SetStats(kFALSE);
  rate->SetMarkerColor(2);
  rate->SetMarkerSize(.5);
  rate->SetMarkerStyle(21);
  rate->GetXaxis()->SetTitle("Time ( OBT in ms )");
  rate->GetXaxis()->CenterTitle();
  rate->GetYaxis()->SetTitle("S4 rate (Hz)");
  rate->GetYaxis()->CenterTitle();
  if(MaxR>500) gPad->SetLogy();
  rate->Draw("9p");

  canvas3->cd(2);
  rateline->SetStats(kFALSE);
  rateline->SetMarkerColor(4);
  rateline->SetMarkerSize(.5);
  rateline->SetMarkerStyle(21);
  rateline->GetXaxis()->SetTitle("Time ( OBT in ms )");
  rateline->GetXaxis()->CenterTitle();
  rateline->GetYaxis()->SetTitle("S4 rate (Hz)");
  rateline->GetYaxis()->CenterTitle();
  rateline->SetMaximum(rate->GetMaximum());
  if(MaxR>500) gPad->SetLogy();
  rateline->Draw("9p");
 

  if (outDir == "./") {
    oss3.str("");
    oss3 << filename.Data() << "_S4_QL_3." << format.Data();
  } else {
    oss3.str("");
    oss3 << outDir.Data() << filename.Data() << "_S4_QL_3." << format.Data();
  }
  canvas3->SaveAs(oss3.str().c_str());
  file->Close();
}


int main(int argc, char* argv[]){
  TString path;
  TString outDir  ="./";
  TString format  ="jpg";
  ULong_t DeltaT  =1;
   
    
 if (argc < 2){
    printf("You have to insert at least the file to analyze \n");
    printf("Try '--help' for more information. \n");
    exit(1);
  }  

  if (!strcmp(argv[1], "--help")){
        printf( "Usage: ND_QL FILE [OPTION] \n");
        printf( "\t --help                 Print this help and exit \n");
        printf( "\t -outDir[path]          Path where to put the output [default ./] \n");
        printf( "\t -format[jpg|gif|ps]    Format for output files [default 'jpg'] \n");
        printf( "\t -DeltaT                Time interval to control time evolution of acquired data in minute [default = 1 min] \n");
        exit(1);
  }

  path=argv[1];

  for (int i = 2; i < argc; i++){
    
    if (!strcmp(argv[i], "-outDir")){
      if (++i >= argc){
        printf( "-outDir needs arguments. \n");
        printf( "Try '--help' for more information. \n");
        exit(1);
      } 
      else{
        outDir = argv[i];
        continue;
      }
    }

    if (!strcmp(argv[i], "-format")){
      if (++i >= argc){
        printf( "-format needs arguments. \n");
        printf( "Try '--help' for more information. \n");
        exit(1);
      }
      else{
        format = argv[i];
        continue;
      }
    }

    if (!strcmp(argv[i], "-DeltaT")){
      if (++i >= argc){
        printf( "-DeltaT needs arguments. \n");
        printf( "Try '--help' for more information. \n");
        exit(1);
      } 
      else{
        DeltaT = atol(argv[i]);
        continue;
      }
    }


  }

   S4_QL(argv[1], outDir, format, DeltaT);

}