examples/test.cpp

// File: test.cpp                                                            //
// Description: example program that implements tests with random particles  //
// This file is part of the SISCone project.                                 //
// For more details, see http://projects.hepforge.org/siscone                //
//                                                                           //
// Copyright (c) 2006 Gavin Salam and Gregory Soyez                          //
//                                                                           //
// 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 2 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, write to the Free Software               //
// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA //
//                                                                           //
// $Revision:: 179                                                          $//
// $Date:: 2007-10-03 13:08:54 -0400 (Wed, 03 Oct 2007)                     $//

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <sys/time.h>
#include <iostream>
#include <math.h>

#include "siscone/ranlux.h"
#include "siscone/momentum.h"
#include "siscone/protocones.h"
#include "siscone/split_merge.h"
#include "siscone/siscone.h"

#define N_default  500
#define R          0.7
#define F          0.75

using namespace std;
using namespace siscone;

int main(int argc, char* argv[]){
  vector<Cmomentum> particles;
  Cmomentum *v;
  double phi=0, eta=0, pt=1;
  unsigned int N;

  unsigned int i, nb_stable;
  FILE *flux;

  if (argc==1){
    cout << "using default number of particles" << endl;
    N = N_default;
  } else {
    sscanf(argv[1], "%d", &N);
    cout << "using " << N << " particles" << endl;
  }    

  // Initialise random number generator
  cout << "initialise random number generator" << endl;
  timeval timestamp;

  gettimeofday(&timestamp, NULL);
  srand(timestamp.tv_usec);


  // build particle list
  cout << "build particle list" << endl;
  flux = fopen("particles.dat", "w+");
  ranlux_init();
  for (i=0;i<N;i++){
    // example 1: 4 particles in lozenge (11 stable cones needed)
    phi = 0.9*sin(2.0*i*M_PI/N);
    if (phi<-M_PI) phi+=2.0*M_PI;
    if (phi> M_PI) phi-=2.0*M_PI;
    eta = 1.1*cos(2.0*i*M_PI/N);

    // example 2: random particles
    eta = -5.0+10.0*rand()/(RAND_MAX+1.0);
    phi = 2.0*M_PI*rand()/(RAND_MAX+1.0);

    // example 3: quadtree needed to check stability
    /*
    switch (i){
    case 0: eta = -0.9; phi = 0.0; break;
    case 1: eta =  0.9; phi = 0.0; break;
    case 2: eta =  0.0; phi = 0.5; break;
    case 3: eta =  0.0; phi = 0.6; break;
    };
    */
    pt = 1.0;

    /*
    switch (i){
    case 0: eta = 0.0; phi = 0.0; pt = 4.0; break;
    case 1: eta = 0.9; phi = 0.0; pt = 1.0; break;
    case 2: eta = 2.3; phi = 0.0; pt = 1.0; break;
    };
    */

    particles.push_back(Cmomentum(pt*cos(phi), pt*sin(phi), pt*sinh(eta), pt*cosh(eta)));

    fprintf(flux, "%le\t%le\n",   particles[i].eta, particles[i].phi);
  }
  fclose(flux);

  cout << "stable cones: initialise engine" << endl;
  Cstable_cones *engine;
  engine = new Cstable_cones(particles);

  // get stable cones
  cout << "stable cones: retreive stable cones" << endl;
  nb_stable = engine->get_stable_cones(R);
  cout << "  found " << nb_stable << " stable cones" << endl;

  // save list of stable cones
  cout << "stable cones: save result" << endl;
  flux = fopen("protocones.dat", "w+");
  for (i=0;i<nb_stable;i++){
    v = &(engine->protocones[i]);
    fprintf(flux, "%le\t%le\n", v->eta, v->phi);
  }
  fclose(flux);

  // splitting and merging
  cout << "split and merge: init" << endl;
  Csplit_merge sm;
  sm.init(particles, &(engine->protocones), R*R);
  delete engine;

  cout << "split and merge: perform" << endl;
  nb_stable = sm.perform(F);
  cout << "  found " << nb_stable << " jets" << endl;
  flux = fopen("jets.dat", "w+");
  sm.save_contents(flux);
  fclose(flux);

  Csiscone *siscone;
  siscone = new Csiscone;
  delete siscone;

  cout << "bye..." << endl;

  return 0;
}

---