Fixing error on gcc 4.5.1
[u/mrichter/AliRoot.git] / TTherminator / Therminator / Hypersurface.cxx
1 #include <cstdlib>
2 #include <cstring>
3 #include "Hypersurface.h"
4
5 Hypersurface::Hypersurface(const char *dirname) {
6 // Read FOHSI.txt file and init parameters
7   char *dirfull;
8   dirfull = (char *) malloc(sizeof(char)*strlen(dirname) + 5);
9   sprintf(dirfull, "%s", dirname);
10   if (strlen(dirfull) > 0) 
11     if (dirfull[strlen(dirfull)-1] != '/')
12       strcat(dirfull, "/");
13   FName = (char *) malloc(sizeof(char)* strlen(dirfull)+ 50);
14   sprintf(FName,"%sFOHSI.txt", dirfull);
15   HSFile = new ifstream;
16   HSFile->open(FName);
17   if(HSFile->is_open()) {
18     HSFile->seekg (0, ios::beg);
19 // phi
20     HSFile->getline(buff,100);
21     Np=atoi(buff);
22     HSFile->getline(buff,100);
23     ip=atof(buff);
24     HSFile->getline(buff,100);
25     fp=atof(buff);
26 // zeta
27     HSFile->getline(buff,100);
28     Nz=atoi(buff);
29     HSFile->getline(buff,100);
30     iz=atof(buff);
31     HSFile->getline(buff,100);
32     fz=atof(buff);
33 // freeze-out temperature
34     HSFile->getline(buff,100);
35     TFO=atof(buff);
36 // hydro initial time
37     HSFile->getline(buff,100);
38     tau0=atof(buff);
39
40     HSFile->close();
41     dp = (fp-ip)/(Np-1);
42     dz = (fz-iz)/(Nz-1);
43
44 // create 2D array
45     aArr   = new double* [Np];
46     vArr   = new double* [Np];
47     dArr   = new double* [Np];
48     DpdArr = new double* [Np];
49     DzdArr = new double* [Np];
50     for(i=0;i<Np;i++) {
51       aArr[i]   = new double [Nz];
52       vArr[i]   = new double [Nz];
53       dArr[i]   = new double [Nz];
54       DpdArr[i] = new double [Nz];
55       DzdArr[i] = new double [Nz];
56     }
57   } else {
58     Np = 0; ip = 0.0; fp = 0.0; dp = 1.0;
59     Nz = 0; iz = 0.0; fz = 0.0; dz = 1.0;
60     TFO = 0.0;
61     aArr   = NULL;
62     vArr   = NULL;
63     dArr   = NULL;
64     DpdArr = NULL;
65     DzdArr = NULL;
66   }
67   delete HSFile;
68
69   // Read FOHSa.txt file and fill array
70   sprintf(FName,"%sFOHSa.txt", dirfull);
71   //  FName = "FOHSa.txt";
72   HSFile = new ifstream;
73   HSFile->open(FName);
74   if(HSFile->is_open()) {
75     HSFile->seekg (0, ios::beg);
76     for(i=0;i<Np;i++) {
77       for(j=0;j<Nz;j++) {
78         HSFile->getline(buff,100);
79         aArr[i][j]=atof(buff);
80       }
81     }
82     HSFile->close();
83   }
84   delete HSFile;
85
86   // Read FOHSv.txt file and fill array
87   sprintf(FName,"%sFOHSv.txt", dirfull);
88   //  FName = "FOHSv.txt";
89   HSFile = new ifstream;
90   HSFile->open(FName);
91   if(HSFile->is_open()) {
92     HSFile->seekg (0, ios::beg);
93     for(i=0;i<Np;i++) {
94       for(j=0;j<Nz;j++) {
95         HSFile->getline(buff,100);
96         vArr[i][j]=atof(buff);
97       }
98     }
99     HSFile->close();
100   }
101   delete HSFile;
102
103   // Read FOHSd.txt file and fill array
104   sprintf(FName,"%sFOHSd.txt", dirfull);
105   //  FName = "FOHSd.txt";
106   HSFile = new ifstream;
107   HSFile->open(FName);
108   if(HSFile->is_open()) {
109     HSFile->seekg (0, ios::beg);
110     for(i=0;i<Np;i++) {
111       for(j=0;j<Nz;j++) {
112         HSFile->getline(buff,100);
113         dArr[i][j]=atof(buff);
114       }
115     }
116     HSFile->close();
117   }
118   delete HSFile;
119
120   // Read FOHSDpd.txt file and fill array
121   sprintf(FName,"%sFOHSDpd.txt", dirfull);
122   //  FName = "FOHSDpd.txt";
123   HSFile = new ifstream;
124   HSFile->open(FName);
125   if(HSFile->is_open()) {
126     HSFile->seekg (0, ios::beg);
127     for(i=0;i<Np;i++) {
128       for(j=0;j<Nz;j++) {
129         HSFile->getline(buff,100);
130         DpdArr[i][j]=atof(buff);
131       }
132     }
133     HSFile->close();
134   }
135   delete HSFile;
136
137   // Read FOHSDzd.txt file and fill array
138   sprintf(FName,"%sFOHSDzd.txt", dirfull);
139   //  FName = "FOHSDzd.txt";
140   HSFile = new ifstream;
141   HSFile->open(FName);
142   if(HSFile->is_open()) {
143     HSFile->seekg (0, ios::beg);
144     for(i=0;i<Np;i++) {
145       for(j=0;j<Nz;j++) {
146         HSFile->getline(buff,100);
147         DzdArr[i][j]=atof(buff);
148       }
149     }
150     HSFile->close();
151   }
152   delete HSFile;
153   free (FName);
154   free (dirfull);
155 }
156
157 Hypersurface::Hypersurface(void) {
158   Hypersurface("./");
159 }
160
161
162 Hypersurface::~Hypersurface(void) {
163   for(i=0;i<Np;i++) {
164     delete[] aArr[i];
165     delete[] vArr[i];
166     delete[] dArr[i];
167     delete[] DpdArr[i];
168     delete[] DzdArr[i];
169   }
170   delete[] aArr;
171   delete[] vArr;
172   delete[] dArr;
173   delete[] DpdArr;
174   delete[] DzdArr;
175 }
176
177 double Hypersurface::fahs(double p, double z) {
178   i=int(p/dp);
179   j=int(z/dz);
180   if(i>=Np-1) i=Np-2;
181   if(j>=Nz-1) j=Nz-2;
182   return (aArr[i][j]   * (i+1-p/dp) + aArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
183          (aArr[i][j+1] * (i+1-p/dp) + aArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
184 }
185
186 double Hypersurface::fvhs(double p, double z) {
187   i=int(p/dp);
188   j=int(z/dz);
189   if(i>=Np-1) i=Np-2;
190   if(j>=Nz-1) j=Nz-2;
191   return (vArr[i][j]   * (i+1-p/dp) + vArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
192          (vArr[i][j+1] * (i+1-p/dp) + vArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
193 }
194
195 double Hypersurface::fdhs(double p, double z) {
196   i=int(p/dp);
197   j=int(z/dz);
198   if(i>=Np-1) i=Np-2;
199   if(j>=Nz-1) j=Nz-2;
200   return (dArr[i][j]   * (i+1-p/dp) + dArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
201          (dArr[i][j+1] * (i+1-p/dp) + dArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
202 }
203
204 double Hypersurface::fDpdhs(double p, double z) {
205   i=int(p/dp);
206   j=int(z/dz);
207   if(i>=Np-1) i=Np-2;
208   if(j>=Nz-1) j=Nz-2;
209   return (DpdArr[i][j]   * (i+1-p/dp) + DpdArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
210          (DpdArr[i][j+1] * (i+1-p/dp) + DpdArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
211 }
212
213 double Hypersurface::fDzdhs(double p, double z) {
214   i=int(p/dp);
215   j=int(z/dz);
216   if(i>=Np-1) i=Np-2;
217   if(j>=Nz-1) j=Nz-2;
218   return (DzdArr[i][j]   * (i+1-p/dp) + DzdArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
219          (DzdArr[i][j+1] * (i+1-p/dp) + DzdArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
220 }