RALICE/icepack/iceconvert/IceRawTWR.cxx

   1 /*******************************************************************************
   2  * Copyright(c) 2003, IceCube Experiment at the South Pole. All rights reserved.
   3  *
   4  * Author: The IceCube RALICE-based Offline 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.
  12  * The authors make no claims about the suitability of this software for
  13  * any purpose. It is provided "as is" without express or implied warranty.
  14  *******************************************************************************/
  15
  16 // $Id$
  17
  18 ///////////////////////////////////////////////////////////////////////////
  19 // Class IceRawTWR
  20 // Conversion of Amanda raw TWR data into IceEvent data structures.
  21 // The code to actually read the TWR raw data structures is an Ralice/IcePack
  22 // implementation of Wolfgang Wagner's (Dortmund University, Germany)
  23 // original read_twr_binary_file.cxx and wf2hit_new.cxx source code.
  24 // The trigger information as encountered in the raw data, is available
  25 // in the IceEvent structure via a device named "Trigger".
  26 // The various triggers (and times) have been stored as different "hits"
  27 // in this "Trigger" device, just like it was done in the IceF2k processor
  28 // for the mu-daq F2K data.
  29 // An indication of the active DAQ system is available in the IceEvent structure
  30 // via a device named "Daq". Here the various daq systems (TWR, Muon, ...)
  31 // from which the actual hits (ADC, LE, TOT) eventually will be composed
  32 // are indicated as "signals" of the device itself.
  33 // This class is derived from AliJob providing a task-based processing
  34 // structure on an event-by-event basis.
  35 // The main object in the job environment is an IceEvent* pointer.
  36 // In case the user has provided sub-tasks, these will be executed
  37 // on an event-by-event basis after the IceEvent structure has been filled
  38 // with the raw TWR data and before the final structures are written out.
  39 // Note that the data structures are only written out if an outputfile has
  40 // been specified via the SetOutputFile memberfunction.
  41 // In case no outputfile has been specified, this class provides a facility
  42 // to investigate/analyse raw TWR data using the Ralice/IcePack analysis tools.
  43 //
  44 // Usage example :
  45 // ---------------
  46 //
  47 // gSystem->Load("ralice");
  48 // gSystem->Load("icepack");
  49 // gSystem->Load("iceconvert");
  50 //
  51 // IceRawTWR q("IceRawTWR","TWR raw data to IcePack data structure conversion");
  52 //
  53 // // Limit the number of entries for testing
  54 // q.SetMaxEvents(10);
  55 //
  56 // // Print frequency to produce a short summary print every printfreq events
  57 // q.SetPrintFreq(1);
  58 //
  59 // // The TWR raw data input filename(s)
  60 // q.AddInputFile("twr_2005_101_009225_0983_57784_57850.dat.twr.to_tape_1");
  61 //
  62 // // Output file for the event structures
  63 // q.SetOutputFile("events.root");
  64 //
  65 // ///////////////////////////////////////////////////////////////////
  66 // // Here the user can specify his/her sub-tasks to be executed
  67 // // on an event-by-event basis after the IceEvent structure
  68 // // has been filled and before the data is written out.
  69 // // Sub-tasks (i.e. a user classes derived from TTask) are entered
  70 // // as follows :
  71 // //
  72 // //    MyXtalk task1("task1","Cross talk correction");
  73 // //    MyClean task2("task2","Hit cleaning");
  74 // //    q.Add(&task1);
  75 // //    q.Add(&task2);
  76 // //
  77 // // The sub-tasks will be executed in the order as they are entered.
  78 // ///////////////////////////////////////////////////////////////////
  79 //
  80 // // Perform the conversion and execute subtasks (if any)
  81 // // on an event-by-event basis
  82 // q.ExecuteJob();
  83 //
  84 //--- Author: Nick van Eijndhoven 12-dec-2006 Utrecht University
  85 //- Modified: NvE $Date$ Utrecht University
  86 ///////////////////////////////////////////////////////////////////////////
  87
  88 #include "IceRawTWR.h"
  89 #include "Riostream.h"
  90
  91 ClassImp(IceRawTWR) // Class implementation to enable ROOT I/O
  92
  93 IceRawTWR::IceRawTWR(const char* name,const char* title) : AliJob(name,title)
  94 {
  95 // Default constructor.
  96 // By default maxevent=-1, split=0, bsize=32000, printfreq=1.
  97
  98  fSplit=0;
  99  fBsize=32000;
 100  fMaxevt=-1;
 101  fPrintfreq=1;
 102  fInfiles=0;
 103  fOutfile=0;
 104 }
 105 ///////////////////////////////////////////////////////////////////////////
 106 IceRawTWR::~IceRawTWR()
 107 {
 108 // Default destructor.
 109
 110  if (fInfiles)
 111  {
 112   delete fInfiles;
 113   fInfiles=0;
 114  }
 115 }
 116 ///////////////////////////////////////////////////////////////////////////
 117 void IceRawTWR::SetMaxEvents(Int_t n)
 118 {
 119 // Set the maximum number of events to be processed.
 120 // n=-1 implies processing of the complete input file, which is the default
 121 // initialisation in the constructor.
 122  fMaxevt=n;
 123 }
 124 ///////////////////////////////////////////////////////////////////////////
 125 void IceRawTWR::SetPrintFreq(Int_t f)
 126 {
 127 // Set the printfrequency to produce info every f events.
 128 // f=1 is the default initialisation in the constructor.
 129  if (f>=0) fPrintfreq=f;
 130 }
 131 ///////////////////////////////////////////////////////////////////////////
 132 void IceRawTWR::SetSplitLevel(Int_t split)
 133 {
 134 // Set the split level for the ROOT data file.
 135 // split=0 is the default initialisation in the constructor.
 136  if (split>=0) fSplit=split;
 137 }
 138 ///////////////////////////////////////////////////////////////////////////
 139 void IceRawTWR::SetBufferSize(Int_t bsize)
 140 {
 141 // Set the buffer size for the ROOT data file.
 142 // bsize=32000 is the default initialisation in the constructor.
 143  if (bsize>=0) fBsize=bsize;
 144 }
 145 ///////////////////////////////////////////////////////////////////////////
 146 void IceRawTWR::AddInputFile(TString name)
 147 {
 148 // Add the name of this TWR raw data input file to the list to be processed.
 149
 150  if (!fInfiles)
 151  {
 152   fInfiles=new TObjArray();
 153   fInfiles->SetOwner();
 154  }
 155
 156  TObjString* s=new TObjString();
 157  s->SetString(name);
 158  fInfiles->Add(s);
 159 }
 160 ///////////////////////////////////////////////////////////////////////////
 161 void IceRawTWR::SetOutputFile(TFile* ofile)
 162 {
 163 // Set the output file for the ROOT data.
 164  if (fOutfile) delete fOutfile;
 165  fOutfile=ofile;
 166 }
 167 ///////////////////////////////////////////////////////////////////////////
 168 void IceRawTWR::SetOutputFile(TString name)
 169 {
 170 // Create the output file for the ROOT data.
 171  if (fOutfile) delete fOutfile;
 172  fOutfile=new TFile(name.Data(),"RECREATE","TWR raw data in IceEvent structure");
 173 }
 174 ///////////////////////////////////////////////////////////////////////////
 175 TFile* IceRawTWR::GetOutputFile()
 176 {
 177 // Provide pointer to the ROOT output file.
 178  return fOutfile;
 179 }
 180 ///////////////////////////////////////////////////////////////////////////
 181 void IceRawTWR::Exec(Option_t* opt)
 182 {
 183 // Job to loop over the specified number of events and convert the
 184 // TWR raw data into the IceEvent structure.
 185 // If maxevents<0 (default) all the entries of the input file
 186 // will be processed.
 187 // Every "printfreq" events a short event summary will be printed.
 188 // The default value is printfreq=1.
 189 // The output will be written on a standard output tree named "T".
 190 //
 191 // Notes :
 192 // -------
 193 // 1) This class is derived from AliJob, allowing a task based processing.
 194 //    After the conversion of a raw data event into an IceEvent structure,
 195 //    the processing of all available sub-tasks (if any) is invoked.
 196 //    This provides an event-by-event (sub)task processing before the
 197 //    final data structures are written out.
 198 // 2) The main object in this job environment is an IceEvent* pointer.
 199
 200  if (!fInfiles)
 201  {
 202   cout << " *IceRawTWR Exec* No data input file(s) specified." << endl;
 203   return;
 204  }
 205
 206  Int_t ninfiles=fInfiles->GetEntries();
 207  if (!ninfiles)
 208  {
 209   cout << " *IceRawTWR Exec* No data input file(s) specified." << endl;
 210   return;
 211  }
 212
 213  TTree* otree=0;
 214  if (fOutfile)
 215  {
 216   otree=new TTree("T","TWR raw data converted to IceEvent structures");
 217   otree->SetDirectory(fOutfile);
 218  }
 219
 220  IceEvent* evt=new IceEvent();
 221  evt->SetTrackCopy(1);
 222  evt->SetDevCopy(1);
 223
 224  // Branch in the tree for the event structure
 225  if (otree) otree->Branch("IceEvent","IceEvent",&evt,fBsize,fSplit);
 226
 227  // Initialise the job working environment
 228  SetMainObject(evt);
 229  if (fOutfile)
 230  {
 231   AddObject(fOutfile);
 232   AddObject(otree);
 233  }
 234
 235  TString inputfile;
 236
 237  cout << " ***" << endl;
 238  cout << " *** Start processing of job " << GetName() << " ***" << endl;
 239  cout << " ***" << endl;
 240  for (Int_t i=0; i<ninfiles; i++)
 241  {
 242   TObjString* sx=(TObjString*)fInfiles->At(i);
 243   if (!sx) continue;
 244   inputfile=sx->GetString();
 245   cout << " TWR raw data input file : " << inputfile.Data() << endl;
 246  }
 247  cout << " Maximum number of events to be processed : " << fMaxevt << endl;
 248  cout << " Print frequency : " << fPrintfreq << endl;
 249  if (fOutfile)
 250  {
 251   cout << " ROOT output file : " << fOutfile->GetName() << endl;
 252   cout << " Output characteristics : splitlevel = " << fSplit << " buffersize = " << fBsize << endl;
 253  }
 254
 255  ListEnvironment();
 256
 257  // Storage of the used parameters in the IceRawTWR device
 258  AliDevice params;
 259  params.SetNameTitle("IceRawTWR","IceRawTWR processor parameters");
 260  params.SetSlotName("Nchannels",1);
 261  params.SetSlotName("Ntriggers",2);
 262  params.SetSlotName("BaselineOffset",3);
 263  params.SetSlotName("NanosecsPerTWRbin",4);
 264  params.SetSignal(float(N_OF_CHANNELS),1);
 265  params.SetSignal(float(N_OF_TRIGGERS),2);
 266  params.SetSignal(float(BASELINE_MEAN_MAGIC),3);
 267  params.SetSignal(float(NSECS_PER_TWR_BIN),4);
 268
 269  // Set DAQ device info
 270  AliDevice daq;
 271  daq.SetName("Daq");
 272  daq.SetSlotName("TWR",1);
 273  daq.SetSignal(1,1);
 274
 275  twr_raw_data_file_t twr_file;
 276  Int_t year,runnum,evtnum;
 277
 278  Int_t error;
 279  UInt_t nhead;
 280
 281  GPS_t gps;
 282  UInt_t gpslow,gpshigh,gpssecs; // The GPS time information
 283  Int_t seconds,nsecs;           // Seconds and nanoseconds since start of the UT year
 284
 285  Int_t nevt=0;
 286  fHeader=0;
 287  for (Int_t ifile=0; ifile<ninfiles; ifile++)
 288  {
 289   TObjString* sx=(TObjString*)fInfiles->At(ifile);
 290   if (!sx) continue;
 291
 292   inputfile=sx->GetString();
 293   if (inputfile=="") continue;
 294
 295   // Open the TWR raw data input file in binary mode
 296   fInput=fopen(inputfile.Data(),"rb");
 297
 298   if (!fInput)
 299   {
 300    cout << " *IceRawTWR Exec* No input file found with name : " << inputfile.Data() << endl;
 301    continue;
 302   }
 303
 304   // Extract info like run number, file number etc... from filename
 305   extract_info_from_filename((char*)inputfile.Data(),&twr_file);
 306
 307   year=twr_file.year;
 308   runnum=twr_file.run_no;
 309
 310   // Initialise the event structure
 311   clear_event(&fEvent);
 312
 313   // Read the file header information
 314   error=read_header_from_file(fInput,&fHeader,&nhead);
 315
 316   if (error || !nhead)
 317   {
 318    cout << " *IceRawTWR Exec* Error in header for input file : " << inputfile.Data() << endl;
 319    continue;
 320   }
 321
 322   // Correct the mapping
 323   update_system(fHeader,runnum);
 324
 325   // Retrieve the actual readout system, threshold and external stop of each OM for these data
 326   fReadout.Set(681);
 327   fThreshold.Set(681);
 328   fExtstop.Set(681);
 329   Int_t ncrates=fHeader->n_crates;
 330   Int_t ntwrs=0;
 331   Int_t omid,readout,thresh,extstop;
 332   for (Int_t icr=0; icr<ncrates; icr++)
 333   {
 334    ntwrs=fHeader->crate[icr]->n_twr;
 335    for (Int_t i_twr=0; i_twr<ntwrs; i_twr++)
 336    {
 337     extstop=fHeader->crate[icr]->twr[i_twr]->ext_stop;
 338     for (Int_t ich=0; ich<CHANNELS_PER_TWR; ich++)
 339     {
 340      omid=fHeader->crate[icr]->twr[i_twr]->om_no[ich];
 341      thresh=fHeader->crate[icr]->twr[i_twr]->threshold[ich];
 342      readout=1+fHeader->crate[icr]->twr[i_twr]->om_is_optical[ich];
 343      if (omid>=0 && omid<=681)
 344      {
 345       fReadout.AddAt(readout,omid-1);
 346       fThreshold.AddAt(thresh,omid-1);
 347       fExtstop.AddAt(extstop,omid-1);
 348      }
 349     }
 350    }
 351   }
 352
 353   while (!read_event(fInput,fHeader,&fEvent))
 354   {
 355    if (fMaxevt>-1 && nevt>=fMaxevt) break;
 356
 357    evtnum=fEvent.eventcounter;
 358
 359    // The GPS telegram info
 360    gps=fEvent.gps;
 361    gpslow=gps.seconds & 0x00FFFFFF; // The low 24 bits of the seconds count
 362    gpshigh=gps.info.bits.seconds;   // The high 8 bits of the seconds count
 363    gpssecs=gpshigh<<24;
 364    gpssecs+=gpslow;
 365
 366    // Seconds and nanoseconds since the start of the UT year
 367    seconds=gpssecs;
 368    nsecs=100*gps.count_10MHz;
 369
 370    // Correction for GPS telegram interpretation in the TWR Daq
 371    if (year<2007) seconds-=24*3600;
 372
 373    // Reset the complete Event structure
 374    evt->Reset();
 375
 376    evt->SetRunNumber(runnum);
 377    evt->SetEventNumber(evtnum);
 378    evt->SetUT(year,0,seconds,nsecs);
 379
 380    evt->AddDevice(params);
 381    evt->AddDevice(daq);
 382
 383    PutTrigger(year);
 384
 385    PutWaveforms(year);
 386
 387    // Invoke all available sub-tasks (if any)
 388    CleanTasks();
 389    ExecuteTasks(opt);
 390
 391    if (fPrintfreq)
 392    {
 393     if (!(nevt%fPrintfreq)) evt->HeaderData();
 394    }
 395
 396    // Write the complete structure to the output Tree
 397    if (otree) otree->Fill();
 398
 399    // Update event counter
 400    nevt++;
 401
 402    // Reset the raw event structure
 403    clear_event(&fEvent);
 404   } // End of event reading loop
 405
 406   // Delete the file header structure
 407   clear_system(fHeader);
 408
 409   if (fMaxevt>-1 && nevt>=fMaxevt) break;
 410
 411  } // End of input file loop
 412
 413  // Flush possible memory resident data to the output file
 414  if (fOutfile) fOutfile->Write();
 415
 416  // Remove the IceEvent object from the environment
 417  // and delete it as well
 418  if (evt)
 419  {
 420   RemoveObject(evt);
 421   delete evt;
 422  }
 423 }
 424 ///////////////////////////////////////////////////////////////////////////
 425 void IceRawTWR::PutWaveforms(Int_t year)
 426 {
 427 // Get the waveform info from the raw data event into the IcePack structure.
 428
 429  IceEvent* evt=(IceEvent*)GetMainObject();
 430  if (!evt) return;
 431
 432  // Get trigger time for TWR time reference
 433  Float_t trigtime=0;
 434  AliDevice* trig=(AliDevice*)evt->GetDevice("Trigger");
 435  if (trig)
 436  {
 437   AliSignal* sx=trig->GetHit("main");
 438   if (sx) trigtime=sx->GetSignal("trig_pulse_le");
 439  }
 440
 441  // Loop over all the waveforms and add the histo(s) to the corresponding OM's
 442  TH1F histo;
 443  Int_t nbins=0;
 444  Float_t xlow=0;
 445  Float_t xup=0;
 446  TString hname;
 447  IceAOM om;
 448  IceAOM* omx=0;
 449  Int_t omid;
 450  Int_t omidmax=680;
 451  Int_t error;
 452  Float_t baseline;
 453  Int_t nfrags;
 454  Int_t firstbin,lastbin;
 455  Int_t extstop;
 456  for (Int_t i=0; i<N_OF_CHANNELS; i++)
 457  {
 458   if (!fEvent.wfm_filled[i]) continue;
 459
 460   omid=i;
 461   if (omid<=0 || omid>omidmax) continue; // Skip trigger channels
 462
 463   // Get corresponding device from the current event structure
 464   omx=(IceAOM*)evt->GetIdDevice(omid);
 465   if (!omx)
 466   {
 467    om.Reset(1);
 468    om.SetUniqueID(omid);
 469    evt->AddDevice(om);
 470    omx=(IceAOM*)evt->GetIdDevice(omid);
 471   }
 472
 473   if (!omx) continue;
 474
 475   extstop=fExtstop.At(omid-1);
 476
 477   // Update readout type and threshold for this OM
 478   omx->AddNamedSlot("READOUT");
 479   omx->SetSignal(float(fReadout.At(omid-1)),"READOUT");
 480   omx->AddNamedSlot("THRESH");
 481   omx->SetSignal(float(fThreshold.At(omid-1)),"THRESH");
 482   omx->AddNamedSlot("EXTSTOP");
 483   omx->SetSignal(float(extstop),"EXTSTOP");
 484
 485   clear_waveform_analysis(&fWform);
 486   error=restore_waveform(fEvent.wfm[i],&fWform,year);
 487
 488   if (error) continue;
 489
 490   nfrags=fWform.n_frag;
 491
 492   for (Int_t ifrag=0; ifrag<nfrags; ifrag++)
 493   {
 494    baseline=fWform.frag_mean[ifrag];
 495
 496    hname="BASELINE-WF";
 497    hname+=omx->GetNwaveforms()+1;
 498    omx->AddNamedSlot(hname);
 499    omx->SetSignal(baseline,hname);
 500
 501    // Fill the waveform histogram with this fragment
 502    hname="OM";
 503    hname+=omid;
 504    hname+="-WF";
 505    hname+=omx->GetNwaveforms()+1;
 506
 507    histo.Reset();
 508    histo.SetName(hname.Data());
 509    nbins=fWform.frag_n_points[ifrag];
 510    firstbin=fWform.frag_begin[ifrag];
 511    lastbin=fWform.frag_end[ifrag];
 512
 513    xlow=fWform.wfm_x[firstbin];
 514    xup=fWform.wfm_x[lastbin];
 515
 516    // Synchronise waveform times with external stop and trigger time
 517    xlow=xlow-float(NSECS_PER_TWR_BIN)*(1024-extstop)+trigtime;
 518    xup=xup-float(NSECS_PER_TWR_BIN)*(1024-extstop)+trigtime;
 519
 520    histo.SetBins(nbins,xlow,xup);
 521
 522    for (Int_t jbin=1; jbin<=nbins; jbin++)
 523    {
 524     histo.SetBinContent(jbin,baseline-fWform.wfm_y[firstbin+jbin-1]);
 525    }
 526
 527    omx->SetWaveform(&histo,omx->GetNwaveforms()+1);
 528   } // End of loop over fragments of this OM
 529  } // End of loop over channels
 530 }
 531 ///////////////////////////////////////////////////////////////////////////
 532 void IceRawTWR::PutTrigger(Int_t year)
 533 {
 534 // Get the trigger info from the raw data event into the IcePack structure.
 535 // Currently only the trigger settings for the years 2005 and 2006 have been
 536 // implemented.
 537 // In addition to the hardware and software triggers as encountered in the
 538 // raw data, an artificial "main" trigger has been introduced.
 539 // This artificial "main" trigger is just an "or" of the standard hard and soft
 540 // triggers (except calibration and random triggers) and serves only to
 541 // provide a generic "main" trigger a la Amanda mu-daq so that the default
 542 // "IceCleanHits" hit cleaning procedure will work correctly.
 543 // The trigger time for the artificial "main" trigger is taken to be the
 544 // time of the earliest hardware trigger pulse. In case there is no hardware
 545 // trigger pulse available, the "main" trigger time is set to 0.
 546 // For other years, only the artificial "main" trigger with a trigger time
 547 // set to 0 will be stored in the IceEvent structure.
 548
 549  // Fill the trigger structure
 550  Int_t error=retrigger(&fEvent,&fTrigger);
 551  if (error) return;
 552
 553  IceEvent* evt=(IceEvent*)GetMainObject();
 554  if (!evt) return;
 555
 556  AliDevice trig;
 557  trig.SetNameTitle("Trigger","Amanda/IceCube event triggers");
 558  AliSignal s;
 559  Float_t trigtime=0;
 560
 561  if (year !=2005 && year != 2006)
 562  {
 563   s.SetName("main");
 564   s.SetTitle("First trigger for TWR time reference");
 565   s.SetUniqueID(0);
 566   s.SetSlotName("trig_pulse_le",1);
 567   s.SetSignal(trigtime,1);
 568   trig.AddHit(s);
 569   // Store the trigger data into the IceEvent structure
 570   evt->AddDevice(trig);
 571   return;
 572  }
 573
 574  // Trigger settings for 2005 and 2006
 575  if (!fTrigger.n_software_trigger && !fTrigger.n_hardware_trigger) return;
 576
 577  TString trignames[N_OF_TRIGGERS]={"m24","m18","string","spase","cal-t0","cal-la","m12",
 578                                    "main-logic","main-or","random","m20-frag","volume"};
 579  Int_t imain=0;
 580  for (Int_t i=0; i<N_OF_TRIGGERS; i++)
 581  {
 582   if (!fTrigger.trigger_active[i]) continue;
 583
 584   s.Reset(1);
 585   s.SetName(trignames[i]);
 586   s.SetUniqueID(i);
 587   trigtime=0;
 588   if (fTrigger.trigger_has_pulse[i]) trigtime=float(fTrigger.trigger_time[i]*NSECS_PER_TWR_BIN);
 589   s.SetSlotName("trig_pulse_le",1);
 590   s.SetSignal(trigtime,1);
 591   trig.AddHit(s);
 592   // Set flag to indicate creation of artificial "main" trigger
 593   if (i!=4 && i!=5 && i!=9) imain=1;
 594  }
 595
 596  // Set the artificial "main" trigger to indicate the first trigger signal
 597  if (imain)
 598  {
 599   s.Reset(1);
 600   s.SetName("main");
 601   s.SetTitle("First trigger for TWR time reference");
 602   s.SetUniqueID(N_OF_TRIGGERS);
 603   s.SetSlotName("trig_pulse_le",1);
 604   trigtime=0;
 605   if (fTrigger.first_trigger>=0) trigtime=float(fTrigger.first_trigger_time*NSECS_PER_TWR_BIN);
 606   s.SetSignal(trigtime,1);
 607   trig.AddHit(s);
 608  }
 609
 610  // Store the trigger data into the IceEvent structure
 611  evt->AddDevice(trig);
 612 }
 613 ///////////////////////////////////////////////////////////////////////////
 614 Int_t IceRawTWR::extract_info_from_filename(char* fname,twr_raw_data_file_t* twr_file)
 615 {
 616   char start_str[20],year_str[20],day_str[20],run_no_str[20],
 617        file_no_str[20],begin_str[20],end_str[20];
 618   char* filename;
 619
 620   filename = strstr(fname, "twr");
 621   if(filename == NULL)
 622   if(strncmp("twr_", start_str, 4))
 623     {
 624       printf("%s\n", filename);
 625       return(ERROR_NOT_VALID_FILENAME);
 626     }
 627
 628   strncpy(start_str, filename, 4);
 629   if(strncmp("twr_", start_str, 4))
 630     {
 631       printf("%s %s\n", filename, start_str);
 632       return(ERROR_NOT_VALID_FILENAME);
 633     }
 634   strncpy(year_str, &filename[4], 4);
 635   twr_file->year = strtol(year_str, 0, 10);
 636
 637   if(twr_file->year==2003)
 638     {
 639       strncpy(day_str, &filename[9], 3);
 640       day_str[3] = '\0';
 641       twr_file->day = strtol(day_str, 0, 10);
 642
 643       strncpy(run_no_str, &filename[13], 4);
 644       run_no_str[4] = '\0';
 645       twr_file->run_no = strtol(run_no_str, 0, 10);
 646
 647       strncpy(file_no_str, &filename[18], 4);
 648       file_no_str[4] = '\0';
 649       twr_file->file_no = strtol(file_no_str, 0, 10);
 650     }
 651
 652   if(twr_file->year==2004)
 653     {
 654       strncpy(day_str, &filename[9], 3);
 655       day_str[3] = '\0';
 656       twr_file->day = strtol(day_str, 0, 10);
 657
 658       strncpy(run_no_str, &filename[13], 4);
 659       run_no_str[4] = '\0';
 660       twr_file->run_no = strtol(run_no_str, 0, 10);
 661
 662       strncpy(file_no_str, &filename[18], 4);
 663       file_no_str[4] = '\0';
 664       twr_file->file_no = strtol(file_no_str, 0, 10);
 665
 666       strncpy(begin_str, &filename[23], 5);
 667       begin_str[5] = '\0';
 668       twr_file->begin = strtol(begin_str, 0, 10);
 669
 670       strncpy(end_str, &filename[29], 5);
 671       end_str[5] = '\0';
 672       twr_file->end = strtol(end_str, 0, 10);
 673     }
 674
 675   if(twr_file->year > 2004)
 676     {
 677       strncpy(day_str, &filename[9], 3);
 678       day_str[3] = '\0';
 679       twr_file->day = strtol(day_str, 0, 10);
 680
 681       strncpy(run_no_str, &filename[13], 6);
 682       run_no_str[6] = '\0';
 683       twr_file->run_no = strtol(run_no_str, 0, 10);
 684
 685       strncpy(file_no_str, &filename[20], 4);
 686       file_no_str[4] = '\0';
 687       twr_file->file_no = strtol(file_no_str, 0, 10);
 688
 689       strncpy(begin_str, &filename[25], 5);
 690       begin_str[5] = '\0';
 691       twr_file->begin = strtol(begin_str, 0, 10);
 692
 693       strncpy(end_str, &filename[31], 5);
 694       end_str[5] = '\0';
 695       twr_file->end = strtol(end_str, 0, 10);
 696     }
 697   return(0);
 698 }
 699 ///////////////////////////////////////////////////////////////////////////
 700 Int_t IceRawTWR::clear_system(sys_config_t* sys)
 701 {
 702 // Deletion of the file header structure.
 703
 704  if (!sys) return 0;
 705
 706  for(Int_t icrate=0; icrate < int(sys->n_crates); icrate++)
 707  {
 708   if (!sys->crate[icrate]) continue;
 709   for(Int_t itwr=0; itwr < int(sys->crate[icrate]->n_twr); itwr++)
 710   {
 711    if (sys->crate[icrate]->twr[itwr]) delete sys->crate[icrate]->twr[itwr];
 712   }
 713   delete sys->crate[icrate];
 714  }
 715  delete sys;
 716  sys=0;
 717  return 0;
 718 }
 719 ///////////////////////////////////////////////////////////////////////////
 720 Int_t IceRawTWR::clear_event(event_t* event_ptr)
 721 {
 722   Int_t i_value;
 723   Int_t *int_ptr = (int*) event_ptr;
 724
 725   for(i_value=0; i_value < int(sizeof(event_t)/sizeof(Int_t)); i_value++)
 726     {
 727       *int_ptr++ = 0;
 728     }
 729   return(0);
 730 }
 731 ///////////////////////////////////////////////////////////////////////////
 732 Int_t IceRawTWR::read_header_from_file(FILE* fin,sys_config_t** system_ptr,UInt_t* header_length)
 733 {
 734   Int_t i_crate, i_twr, i_channel;
 735   UInt_t count_twr_in_system = 0;
 736   UInt_t dummy;
 737
 738   sys_config_t *sys;
 739
 740   // allocating memory for sys_config structure
 741   sys = (sys_config_t*) malloc( sizeof(sys_config_t) );
 742
 743   fread(&dummy,sizeof(UInt_t),1,fin); // Header Begin Mark
 744
 745   fread(header_length,sizeof(UInt_t),1,fin);   // Length of header
 746   fread(&sys->clockdiv,sizeof(UInt_t),1,fin);
 747   fread(&sys->n_crates,sizeof(UInt_t),1,fin);
 748
 749   if( (sys->n_crates > MAX_N_CRATES) || (sys->n_crates < 0) )
 750     return(ERROR_TOO_MANY_CRATES);
 751
 752   for(i_crate=0; i_crate < int(sys->n_crates); i_crate++)
 753     {
 754       sys->crate[i_crate] =
 755         (crate_config_t*) malloc( sizeof(crate_config_t) );
 756
 757       fread(&sys->crate[i_crate]->vme_base_bridge,sizeof(UInt_t),1,fin);
 758       fread(&sys->crate[i_crate]->vme_base_100MHz,sizeof(UInt_t),1,fin);
 759       fread(&sys->crate[i_crate]->base_gps,sizeof(UInt_t),1,fin);
 760       fread(&sys->crate[i_crate]->n_twr,sizeof(UInt_t),1,fin);
 761
 762       if( (sys->crate[i_crate]->n_twr > MAX_N_TWR_PER_CRATE)
 763           || (sys->crate[i_crate]->n_twr < 0) )
 764         return(ERROR_TOO_MANY_TWRS);
 765
 766       for(i_twr=0; i_twr < int(sys->crate[i_crate]->n_twr); i_twr++)
 767         {
 768           sys->crate[i_crate]->twr[i_twr] =
 769             (twr_config_t*) malloc( sizeof(twr_config_t) );
 770           count_twr_in_system++;
 771           fread(&sys->crate[i_crate]->twr[i_twr]->base,
 772                sizeof(UInt_t),1,fin);
 773           fread(&sys->crate[i_crate]->twr[i_twr]->id,
 774                sizeof(UInt_t),1,fin);
 775
 776           sys->crate[i_crate]->twr[i_twr]->id
 777             = sys->crate[i_crate]->twr[i_twr]->id - 0x10; /* Correct */
 778
 779
 780           fread(&dummy,sizeof(UInt_t),1,fin); /* stat_reg */
 781           fread(&sys->crate[i_crate]->twr[i_twr]->mod_id,
 782                sizeof(UInt_t),1,fin);
 783           fread(&dummy,sizeof(UInt_t),1,fin); /* acq_ctrl */
 784           fread(&sys->crate[i_crate]->twr[i_twr]->ext_start,
 785                sizeof(UInt_t),1,fin);
 786           fread(&sys->crate[i_crate]->twr[i_twr]->ext_stop,
 787                sizeof(UInt_t),1,fin);
 788           fread(&dummy,sizeof(UInt_t),1,fin); /* evtconfig */
 789
 790           for(i_channel = 0; i_channel < CHANNELS_PER_TWR; i_channel++)
 791             {
 792               fread(&sys->crate[i_crate]->twr[i_twr]->om_no[i_channel],
 793                    sizeof(UInt_t),1,fin);
 794             }
 795
 796           for(i_channel = 0; i_channel < CHANNELS_PER_TWR; i_channel++)
 797             {
 798               fread(&sys->crate[i_crate]->twr[i_twr]->om_is_optical[i_channel],
 799                    sizeof(UInt_t),1,fin);
 800             }
 801
 802           for(i_channel = 0; i_channel < CHANNELS_PER_TWR; i_channel++)
 803             {
 804               fread(&sys->crate[i_crate]->twr[i_twr]->baseline[i_channel],
 805                    sizeof(UInt_t),1,fin);
 806             }
 807
 808           for(i_channel = 0; i_channel < CHANNELS_PER_TWR; i_channel++)
 809             {
 810               fread(&sys->crate[i_crate]->twr[i_twr]->threshold[i_channel],
 811                    sizeof(UInt_t),1,fin);
 812             }
 813
 814           sys->twr_field[(i_crate * 0x10) + i_twr]
 815               = sys->crate[i_crate]->twr[i_twr];
 816
 817           /* Bug fix needed */
 818           for(i_channel=0; i_channel < 8; i_channel++)
 819           {
 820              if( sys->crate[i_crate]->twr[i_twr]->om_no[i_channel] == 9000 )
 821                  sys->crate[i_crate]->twr[i_twr]->om_no[i_channel]
 822                    = N_OF_CHANNELS - 1;
 823           }
 824         }
 825     }
 826
 827   // Set number of TWRs in system
 828   sys->n_twr = count_twr_in_system;
 829
 830   *system_ptr = sys;
 831   return(0);
 832 }
 833 ///////////////////////////////////////////////////////////////////////////
 834 Int_t IceRawTWR::update_system(sys_config_t* sys,Int_t run_number)
 835 {
 836   Int_t i_crate, i_twr, i_channel;
 837
 838   // Data for bug fix 1 by Andreas
 839   UInt_t om_no_r1[CHANNELS_PER_TWR]
 840     = {111, 112, 113, 114, 115, 116, 39, 118};
 841   UInt_t om_is_optical_r1[CHANNELS_PER_TWR]
 842      = {0, 0, 0, 0, 0, 0, 0, 0};
 843   UInt_t threshold_r1[CHANNELS_PER_TWR]
 844     = {50, 50, 50, 50, 50, 50, 80, 50};
 845
 846   // Data for bugfix 2 by Timo
 847   UInt_t om_no_r2[CHANNELS_PER_TWR]
 848     = {473, 484, 485, 486, 487, 475, 490, 491};
 849   UInt_t om_is_optical_r2[CHANNELS_PER_TWR]
 850      = {1, 1, 1, 1, 1, 1, 1, 1};
 851   UInt_t threshold_r2[CHANNELS_PER_TWR]
 852     = {15, 50, 55, 40, 15, 23, 15, 15};
 853
 854   // Data for bugfix 3 by Wolfgang
 855   // Old (=incorrect) OM config : {183, 184, 185, 38, 187, 188, 189, 190}
 856   UInt_t om_no_r3[CHANNELS_PER_TWR]={345, 184, 185, 38, 187, 188, 189, 190};
 857   UInt_t om_is_optical_r3[CHANNELS_PER_TWR]={1, 0, 0, 0, 0, 0, 0, 0};
 858   UInt_t threshold_r3[CHANNELS_PER_TWR]={20, 50, 50, 80, 50, 50, 50, 50};
 859
 860   // Old (=incorrect) OM config : {345, 346, 454, 450, 635, 10000, 10000, 10000}
 861   UInt_t om_no_r4[CHANNELS_PER_TWR]={183, 346, 454, 450, 635, 10000, 10000, 10000};
 862   UInt_t om_is_optical_r4[CHANNELS_PER_TWR]={0, 1, 1, 1, 1, 1, 1, 1};
 863   UInt_t threshold_r4[CHANNELS_PER_TWR]={50, 20, 20, 20, 20, 500, 500, 500};
 864
 865   // Bugfix 1 Andreas Bug
 866
 867   //  By accident this TWR was counted twice in TWR.cnf
 868   //  as Crate 0 TWR 7 and Crate 4 TWR 7
 869   //  from run 9153 up to run 9841 (incl.)
 870   //  TWR_OM          639     642     1       9       10      11      12      30
 871   //  OPTICAL         0       0       0       0       0       0       0       0
 872   //  TWR_BASELINE    110     120     110     140     150     160     170     180
 873   //  TWR_THRESHOLD   50      50      80      80      80      80      80      80
 874
 875   //  Crate 4 TWR 7 should be replaced with this TWR
 876   //  TWR_OM          111     112     113     114     115     116     39      118
 877   //  OPTICAL         0       0       0       0       0       0       0       0
 878   //  TWR_BASELINE    110     120     130     140     150     160     170     180
 879   //  TWR_THRESHOLD   50      50      50      50      50      50      80      50
 880
 881   // Begin season 2005 13-feb-2005
 882   // Timo corrected TWR.cnf on 05-apr-2006 after run 9841
 883   if (run_number>=9153 && run_number<=9841)
 884   {
 885    i_crate = 4;
 886    i_twr = 7;
 887    for(i_channel = 0; i_channel < CHANNELS_PER_TWR; i_channel++)
 888    {
 889     sys->crate[i_crate]->twr[i_twr]->om_no[i_channel]=om_no_r1[i_channel];
 890     sys->crate[i_crate]->twr[i_twr]->om_is_optical[i_channel]=om_is_optical_r1[i_channel];
 891     sys->crate[i_crate]->twr[i_twr]->threshold[i_channel]=threshold_r1[i_channel];
 892    }
 893   }
 894
 895   // Bugfix 2 Timos Bug
 896
 897   //  By accident this TWR was counted twice in TWR.cnf
 898   //  as Crate 0 TWR 1 and Crate 5 TWR b
 899   //  from run 9153 up to run 9188 (incl.)
 900
 901   //  TWR_OM          492     493     495     496     497     499     500     501
 902   //  OPTICAL         1       1       1       1       1       1       1       1
 903   //  TWR_BASELINE    110     120     130     140     150     160     170     180
 904   //  TWR_THRESHOLD   16      45      25      42      35      46      15      15
 905
 906   //  Crate 5 TWR b should be corrected to
 907   //  TWR_OM          473     484     485     486     487     475     490     491
 908   //  OPTICAL         1       1       1       1       1       1       1       1
 909   //  TWR_BASELINE    4000    120     130     140     150     4000    170     180
 910   //  TWR_THRESHOLD   15      50      55      40      15      23      15      15
 911
 912  // Begin season 2005 : 13-feb-2005
 913  // Timo corrected TWR.cnf on 15-mar-2005 (= day 74) after run 9188
 914   if (run_number>=9153 && run_number<=9188)
 915   {
 916    i_crate = 5;
 917    i_twr = 0xb;
 918    for(i_channel=0; i_channel<CHANNELS_PER_TWR; i_channel++)
 919    {
 920     sys->crate[i_crate]->twr[i_twr]->om_no[i_channel]=om_no_r2[i_channel];
 921     sys->crate[i_crate]->twr[i_twr]->om_is_optical[i_channel]=om_is_optical_r2[i_channel];
 922     sys->crate[i_crate]->twr[i_twr]->threshold[i_channel]=threshold_r2[i_channel];
 923    }
 924   }
 925
 926  // Bugfix 3 by Wolfgang Wagner : Mismatch of OM 345 and 183
 927  // Begin season 2005 : 13-feb-2005
 928  // Wolfgang corrected TWR.cnf on 15-aug-2005 after run 9988
 929  if (run_number>=9153 && run_number<=9988)
 930  {
 931   i_crate=5;
 932   i_twr=4;
 933   for(i_channel=0; i_channel<CHANNELS_PER_TWR; i_channel++)
 934   {
 935    sys->crate[i_crate]->twr[i_twr]->om_no[i_channel]=om_no_r3[i_channel];
 936    sys->crate[i_crate]->twr[i_twr]->om_is_optical[i_channel]=om_is_optical_r3[i_channel];
 937    sys->crate[i_crate]->twr[i_twr]->threshold[i_channel]=threshold_r3[i_channel];
 938   }
 939
 940   i_twr=5;
 941   for(i_channel=0; i_channel<CHANNELS_PER_TWR; i_channel++)
 942   {
 943    sys->crate[i_crate]->twr[i_twr]->om_no[i_channel]=om_no_r4[i_channel];
 944    sys->crate[i_crate]->twr[i_twr]->om_is_optical[i_channel]=om_is_optical_r4[i_channel];
 945    sys->crate[i_crate]->twr[i_twr]->threshold[i_channel]=threshold_r4[i_channel];
 946   }
 947  }
 948
 949  return(0);
 950 }
 951 ///////////////////////////////////////////////////////////////////////////
 952 Int_t IceRawTWR::read_event(FILE* fin,sys_config_t* sys,event_t* event_ptr)
 953 {
 954     Int_t i_wfm;
 955     UInt_t length_of_event_block;
 956
 957     Int_t n_twr, n_of_waveforms_in_event, read_number;
 958     UInt_t length_wfm[CHANNELS_PER_TWR];
 959     UInt_t dummy, channel_no, om_no, twr_no;
 960
 961     // Reset waveform filled register
 962     memset(&event_ptr->wfm_filled[0], 0, sizeof(UInt_t) * N_OF_CHANNELS);
 963
 964     if( !fread(&dummy,sizeof(UInt_t),1,fin) ) return(1);
 965
 966     if(dummy != 0xbbbbbbbb)
 967     {
 968         printf("Wrong event begin mark %x\n", dummy);
 969         while( (dummy !=0xbbbbbbbb)
 970                && (fread(&dummy,sizeof(UInt_t),1,fin) != 0) )
 971           {;//printf("dummy:%x\n", dummy);
 972           }
 973     }
 974     if( !fread(&length_of_event_block,sizeof(UInt_t),1,fin) ) return(1);
 975     if( !fread(&event_ptr->eventcounter,sizeof(UInt_t),1,fin) ) return(1);
 976     if( !fread(&event_ptr->which_trigger,sizeof(UInt_t),1,fin) ) return(1);
 977     if( !fread(&event_ptr->gps,sizeof(GPS_t),1,fin) ) return(1);
 978
 979     // --reading waveforms from TWR blocks
 980     n_twr = 0;
 981     while(n_twr < int(sys->n_twr))
 982     {
 983         // --read TWR header
 984         if( !fread(&dummy,sizeof(UInt_t),1,fin) ) return(1);
 985         if(dummy != 0xffffffff)
 986         {printf("Wrong twr begin mark %x\n", dummy); return(2);}
 987         if( !fread(&twr_no,sizeof(UInt_t),1,fin) ) return(1);
 988
 989         // nur voruebergehend !!
 990         twr_no -= 0x10;
 991
 992         if( !fread(&event_ptr->twr[twr_no].timestamp,sizeof(UInt_t),1,fin) )
 993             return(1);
 994         if( !fread(&n_of_waveforms_in_event,sizeof(UInt_t),1,fin) )
 995             return(1);
 996         event_ptr->twr[twr_no].n_wfm = n_of_waveforms_in_event;
 997
 998         for(i_wfm=0; i_wfm < n_of_waveforms_in_event; i_wfm++)
 999         {
1000             if( !fread(&length_wfm[i_wfm],sizeof(UInt_t),1,fin) ) return(1);
1001         }
1002
1003         // read waveforms
1004         for(i_wfm=0; i_wfm < n_of_waveforms_in_event; i_wfm++)
1005         {
1006             if(length_wfm[i_wfm] != 0)
1007             {
1008               if( !fread(&channel_no,sizeof(UInt_t),1,fin) ) return(1);
1009               if(sys->twr_field[twr_no]->om_no[channel_no]
1010                    < N_OF_CHANNELS)
1011                     om_no = sys->twr_field[twr_no]->om_no[channel_no];
1012                 else
1013                     om_no = N_OF_CHANNELS-1;
1014
1015                 /* Fix needed */
1016
1017                 event_ptr->twr_id_of_om[om_no] = twr_no;
1018
1019                 read_number = fread(&event_ptr->wfm[om_no],
1020                                    length_wfm[i_wfm]-sizeof(UInt_t),1,fin);
1021                 event_ptr->wfm_filled[om_no] = 1;
1022                 if( !read_number ) return(1);
1023
1024                 // read_number correction for usage of fread() instead of read()
1025                 read_number*=length_wfm[i_wfm]-sizeof(UInt_t);
1026
1027                 if( read_number != int(length_wfm[i_wfm]-sizeof(UInt_t)) )
1028                 {
1029                   cout << " read_number : " << read_number
1030                        << " length_wfm["<<i_wfm<<"] : " << length_wfm[i_wfm]
1031                        << " sizeof(UInt_t) : " << sizeof(UInt_t) << endl;
1032                   return(2);
1033                 }
1034             }
1035         }
1036         n_twr++;
1037     } // end while n_twr
1038     return(0);
1039 }
1040 ///////////////////////////////////////////////////////////////////////////
1041 Int_t IceRawTWR::retrigger(event_t* ev,trigger_hits_t* trig)
1042 {
1043 // Returns the active trigger(s)
1044
1045  // Initialise the trigger_hits_t structure with zeroes
1046  memset(trig, 0, sizeof(trigger_hits_t) );
1047
1048  // Obtain the software trigger info
1049  trig->n_software_trigger=0;
1050  for(Int_t itrigger=0; itrigger<N_OF_TRIGGERS; itrigger++)
1051  {
1052   if(ev->which_trigger & trigger_bits[itrigger])
1053   {
1054    //printf("SetTrigger %i\n", i_trigger);
1055    trig->trigger_active[itrigger]=1;
1056    trig->n_software_trigger++;
1057   }
1058   else
1059   {
1060    trig->trigger_active[itrigger]=0;
1061   }
1062  }
1063
1064  // Obtain the hardware trigger info
1065  trig->n_hardware_trigger=0;
1066  trig->first_trigger_time=10000000;
1067  trig->first_trigger=-1;
1068
1069  for(Int_t jtrigger=0; jtrigger<N_OF_TRIGGERS; jtrigger++)
1070  {
1071   if(!trigger_channel[jtrigger]) continue;
1072
1073   if(ev->wfm_filled[trigger_channel[jtrigger]])
1074   {
1075    trig->trigger_active[jtrigger]=1;
1076    trig->trigger_time[jtrigger]=(ev->wfm[trigger_channel[jtrigger]].value[2] & 0xfff);
1077    trig->trigger_has_pulse[jtrigger]=1;
1078    if (trig->trigger_time[jtrigger] < trig->first_trigger_time)
1079    {
1080     trig->first_trigger_time=trig->trigger_time[jtrigger];
1081     trig->first_trigger=jtrigger;
1082    }
1083    trig->n_hardware_trigger++;
1084   }
1085  }
1086  return 0;
1087 }
1088 ///////////////////////////////////////////////////////////////////////////
1089 Int_t IceRawTWR::clear_waveform_analysis(waveform_analyse_t* wfm_om)
1090 {
1091   Int_t i_value, i_frag, i_edge, i_peak;
1092
1093   if(wfm_om == 0) return(1);
1094
1095   // output from analysis
1096   wfm_om->n_frag = 0;
1097   for(i_frag=0; i_frag < MAX_N_OF_FRAGS; i_frag++)
1098     {
1099       wfm_om->frag_n_points[i_frag] = 0;
1100       wfm_om->frag_begin[i_frag] = 0;
1101       wfm_om->frag_end[i_frag] = 0;
1102       wfm_om->frag_mean[i_frag] = 0;
1103       wfm_om->frag_begin_time[i_frag] = 0;
1104     }
1105
1106   wfm_om->n_peak = 0;
1107   for(i_peak=0; i_peak < MAX_N_OF_PEAKS; i_peak++)
1108     {
1109       wfm_om->peak_begin[i_peak] = 0;
1110       wfm_om->peak_end[i_peak] = 0;
1111       wfm_om->peak_max[i_peak] = 0;
1112       wfm_om->peak_TDC_edge[i_peak] = 0;
1113       wfm_om->peak_local_minimum[i_peak] = 0;
1114       wfm_om->crosstalk_charge_n_value[i_peak] = 0;
1115       wfm_om->peak_in_fragment[i_peak] = 0;
1116
1117       wfm_om->peak_mean[i_peak] = 0.0;
1118
1119       wfm_om->peak_m[i_peak] = 0.0;
1120       wfm_om->peak_b[i_peak] = 0.0;
1121       wfm_om->peak_t0[i_peak] = 0.0;
1122       wfm_om->peak_begin_time[i_peak] = 0.0;
1123       wfm_om->peak_charge[i_peak] = 0.0;
1124       wfm_om->peak_height[i_peak] = 0.0;
1125       wfm_om->fitted_amplitude[i_peak] = 0.0;
1126       wfm_om->fitted_TOT[i_peak] = 0.0;
1127       wfm_om->crosstalk_charge[i_peak] = 0.0;
1128       wfm_om->crosstalk_slope[i_peak] = 0.0;
1129     }
1130
1131   wfm_om->n_point = 0;
1132   wfm_om->wfm_min = 4095;
1133   wfm_om->wfm_max = 0;
1134   wfm_om->b_out_of_range = 0;
1135
1136   for(i_value=0; i_value < 1024; i_value++)
1137     {
1138       wfm_om->wfm_x[i_value] = 0;
1139       wfm_om->wfm_y[i_value] = 0;
1140     }
1141
1142   wfm_om->n_tdc_edges = 0;
1143   for(i_edge=0; i_edge < MAX_N_OF_TDC_EDGES; i_edge++)
1144     {
1145       wfm_om->leading_edge[i_edge] = 0.0;
1146       wfm_om->falling_edge[i_edge] = 0.0;
1147       wfm_om->identified_twr_hit[i_edge] = -1;
1148     }
1149
1150  return(0);
1151 }
1152 ///////////////////////////////////////////////////////////////////////////
1153 Int_t IceRawTWR::restore_waveform(waveform_t f_wfm,waveform_analyse_t* wfm_om,Int_t year)
1154 {
1155     UShort_t wfm_length, mean;
1156     static UShort_t tmp_wf[2000];
1157
1158     Int_t debug = 0;
1159     Int_t fragment_start = 0;
1160     Int_t frag_count = 0;     // position in current fragment
1161     Int_t n_position = 0;     // position in displayed waveform
1162     UInt_t  n_word = 2;      // position in featured waveform
1163     Int_t n_fragment = 0;     // actual fragment
1164     Int_t b_wrong_value = 0;
1165
1166     UShort_t assumed_frag_begin, last_value; /* bug in eventbuilder */
1167
1168     wfm_om->wfm_min = 4095.0;
1169     wfm_om->wfm_max = 0.0;
1170
1171     if( (f_wfm.value[0] & 0xf000) != 0xf000 ) return(1);
1172     wfm_length = (f_wfm.value[0] & 0xfff)/2;
1173
1174     mean = f_wfm.value[1] + BASELINE_MEAN_MAGIC;
1175     while( ((f_wfm.value[n_word] & 0xf000) == 0x4000) &&
1176            (n_word < wfm_length) &&
1177            (n_fragment < MAX_N_OF_FRAGS) )
1178     {
1179       fragment_start = f_wfm.value[n_word] & 0xfff;
1180       n_word++;
1181       wfm_om->frag_begin_time[n_fragment]
1182         = fragment_start * NSECS_PER_TWR_BIN;
1183       wfm_om->frag_begin[n_fragment] = n_position;
1184       wfm_om->frag_mean[n_fragment] = mean;
1185
1186       b_wrong_value = 0;
1187       frag_count = 0;
1188
1189       while( ((f_wfm.value[n_word] & 0xf000) != 0x2000) &&
1190              ((f_wfm.value[n_word] & 0xf000) != 0x4000) &&/*Reconstructable*/
1191              !b_wrong_value && /* Buggy */
1192              (n_word < wfm_length) )
1193         {
1194           if(year > 2004)
1195             {
1196               /* 2005 2006 data */
1197               if(frag_count == 0)
1198                 {
1199                   tmp_wf[n_word] = f_wfm.value[n_word] + mean;
1200                   wfm_om->wfm_y[n_position] = (float) tmp_wf[n_word];
1201                   wfm_om->wfm_x[n_position] = (float)
1202                     wfm_om->frag_begin_time[n_fragment]
1203                     + (frag_count * NSECS_PER_TWR_BIN);
1204                 }
1205               else if(frag_count == 1)
1206                 {
1207                   tmp_wf[n_word] = f_wfm.value[n_word] + tmp_wf[n_word-1];
1208                   wfm_om->wfm_y[n_position] = (float) tmp_wf[n_word];
1209                   wfm_om->wfm_x[n_position] = (float)
1210                     wfm_om->frag_begin_time[n_fragment]
1211                     + (frag_count * NSECS_PER_TWR_BIN);
1212                 }
1213               else
1214                 {
1215                   tmp_wf[n_word] =
1216                     2*tmp_wf[n_word-1] + f_wfm.value[n_word];
1217                   tmp_wf[n_word] -= tmp_wf[n_word-2];
1218
1219                   wfm_om->wfm_y[n_position] = (float) tmp_wf[n_word];
1220                   wfm_om->wfm_x[n_position] = (float)
1221                     wfm_om->frag_begin_time[n_fragment]
1222                     + (frag_count * NSECS_PER_TWR_BIN);
1223                 }
1224
1225
1226               /*
1227                 Hack for wrongly merged overlapping fragments
1228               */
1229               if(tmp_wf[n_word] > 0x1fff)
1230                 {
1231                   /* BUG FIXXXX                                         */
1232                   /* assume that fragment merge in eventbuilder caused  */
1233                   /* problem two fragments overlap in EXACTLY ONE point */
1234                   /* and are merged first point of the added part of    */
1235                   /* the fragment is encoded using the former fragment  */
1236                   /* start as a data point                              */
1237
1238                   last_value         = tmp_wf[n_word-1];
1239                   assumed_frag_begin = 0x4000 + fragment_start + frag_count;
1240                   tmp_wf[n_word] =  f_wfm.value[n_word] + 2 * last_value;
1241                   tmp_wf[n_word] -= assumed_frag_begin;
1242                   wfm_om->wfm_y[n_position] = (float) tmp_wf[n_word];
1243
1244                   /* Look if value is still buggy */
1245                   if(tmp_wf[n_word] > 0x1fff) b_wrong_value = 1;
1246
1247                   debug = ERROR_MISS_FRAG_STOP;
1248                 }
1249             } /* end year >= 2005 */
1250           else
1251             {
1252               /* 2003 2004 data */
1253               wfm_om->wfm_y[n_position] = (float) f_wfm.value[n_word];
1254               wfm_om->wfm_x[n_position] = (float)
1255                 wfm_om->frag_begin_time[n_fragment]
1256                 + (frag_count * NSECS_PER_TWR_BIN);
1257             } /* end year 2003 2004 */
1258
1259           /* Set min and max Y */
1260
1261           if(wfm_om->wfm_y[n_position] > wfm_om->wfm_max)
1262             wfm_om->wfm_max = wfm_om->wfm_y[n_position];
1263           if(wfm_om->wfm_y[n_position] < wfm_om->wfm_min)
1264             wfm_om->wfm_min = wfm_om->wfm_y[n_position];
1265
1266           n_position++;
1267           n_word++;
1268           frag_count++;
1269         }
1270
1271         if((f_wfm.value[n_word] & 0xf000) == 0x2000) /* Normal wavf */
1272           {
1273
1274             wfm_om->frag_end[n_fragment] = n_position - 1;
1275             wfm_om->frag_n_points[n_fragment] =
1276               wfm_om->frag_end[n_fragment]
1277               - wfm_om->frag_begin[n_fragment] + 1;
1278             wfm_om->n_point += wfm_om->frag_n_points[n_fragment];
1279             n_word++;
1280           }
1281         else
1282           return(ERROR_CORRUPTED_WF);
1283
1284         n_fragment++;
1285     } /* end while fragment */
1286
1287
1288     wfm_om->n_frag = n_fragment;
1289     if( !(n_word & 0x1) ) n_word++;
1290
1291     if(n_fragment >= MAX_N_OF_FRAGS) return(ERROR_MAX_N_FRAGMENTS_EXCEEDED);
1292
1293
1294     // Hack to get rid of last value of waveform always set to 0
1295     if (wfm_om->wfm_y[wfm_om->n_point] == 0.0)
1296     {
1297      // erase last point of waveform
1298      wfm_om->n_point--;
1299
1300      // Shorten last pulse if necessary
1301      // if( wfm_om.peak_end[wfm_om.n_peak-1]
1302      //     == wfm_om.frag_end[wfm_om.n_frag-1] )
1303      //   wfm_om.peak_end[wfm_om.n_peak-1]--;
1304
1305      // Shorten last fragment
1306      wfm_om->frag_n_points[wfm_om->n_frag-1]--;
1307      wfm_om->frag_end[wfm_om->n_frag-1]--;
1308
1309      wfm_om->wfm_min = 4095.0;
1310      wfm_om->wfm_max = 0.0;
1311      for (Int_t i_value=0; i_value < wfm_om->n_point; i_value++)
1312      {
1313       if (wfm_om->wfm_y[i_value] > wfm_om->wfm_max) wfm_om->wfm_max=wfm_om->wfm_y[i_value];
1314       if (wfm_om->wfm_y[i_value] < wfm_om->wfm_min) wfm_om->wfm_min=wfm_om->wfm_y[i_value];
1315      }
1316     }
1317
1318  return(debug);
1319 }
1320 ///////////////////////////////////////////////////////////////////////////