Projects/Clusterer.C at master · skyriacoCMS/Projects · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
#include <iostream>
#include <fstream>
#include <math.h>
#include <cstdlib>
#include <vector>

using namespace std;

void Cluster(int x_i, int y_j, vector<double> posx,vector<double> posy, double dens1, vector<int> &x_1 , vector<int> &y_1 ){


  int Grid2[10][10] = {};

  for(int i  = 0; i < posx.size(); i ++){

    double x = posx[i];
    double y = posy[i];

    if(int(x/100) != x_i)continue;
    if(int(y/100) != y_j)continue;

    int xi = int(int(x) %100)/10;
    int yi = int(int(y) %100)/10;


    Grid2[xi][yi] ++;
  }

  for(int i = 0; i < 10 ; i++){
    for(int j = 0; j <10 ; j++){


      double dens2 = Grid2[i][j]/100.;

      if(dens1 <  dens2){
	cout<<"found cluster in sub area:"<<x_i<<"  "<<y_j<<endl;
	x_1.push_back(i);
	y_1.push_back(j);
      }


      }
  }

}


int main(){


  cout<<"--------------Cluster Finder------------"<<endl;

  ofstream out2;

  out2.open("Haloes.dat");


  ifstream in;
  in.open("Cluster.dat",ios::in);

  double number;
  vector<double> posx,posy;

  int i = 0;

  //feeding in the data
  while(!in.eof()){

    in>>number;

    //cout<<number<<endl;


    if( i%2 == 0 )posx.push_back(number);
    if( i%2 == 1) posy.push_back(number);
    //cout<<i%2<<endl;

    i++;
  }

  //general density of space


  double dens0 = posx.size()/(1000.*1000.);


  int Grid [10][10] = {};


  for (int i = 0; i  < posx.size(); i++){

    double x = posx[i];
    double y = posy[i];

    int ix = int(x/100);
    int iy = int(y/100);

    Grid[ix][iy] ++;

  }


  for(int i  = 0; i < 10; i++){
    for(int j = 0; j < 10 ; j++){
      vector<int> x_1;
      vector<int> y_1;

      double dens1 = Grid[i][j]*100./(1000000);
      if(dens1 > dens0) Cluster(i,j, posx,posy,dens1,x_1,y_1 );
      for(int ik = 0; ik < x_1.size();ik++){

	  out2<< i*100 + x_1[ik]*10 <<"  "<< j*100 + y_1[ik]*10<<endl;

      }

    }
  }

  return 0;


  posx.resize(0);
  posy.resize(0);

}