Allink  v0.1
UsaMatematica.cpp
1 #include "../../include/Matematica.h"
2 #include "../../include/VarDatFile.h"
3 
4 #include <iostream>
5 using namespace std;
6 
7 int main(int argc,char **argv){
8  FILE *SEGNALE;
9  char Uscita[256];
10  char Comando[50];
11  char Executex[60];
12  char Executey[60];
13  double *Dis = new double;
14  int CoordX = 0;
15  int NBin = 200;
16  int NMax=0;
17  int NVar =0;
18  int CoordY = 1;
19  int CoordY1 = 1;
20  int CoordX1 = 1;
21  int NMedia=10;
22  int LogLog = 0;
23  int IfNorm = 0;
24  int IfBorder = 0;
25  int NSample = 10;
26  double *Intervalli;
27  double *dInter;
28  //vector <double> st;
29  double *sp;
30  double *sw;
31  double Border[4] = {0.,0.,0.,0.};
32  double Param = 1.;
33  VarDatFile *v1;
34  Matematica *m1 = new Matematica();
35  MOMENTI Mom;
36  sprintf(Uscita,"Risposta.dat");
37  // for(int i=0;i<argc;i++){
38  int NFile = 0;
39  int *ArgFile = (int *)calloc(argc,sizeof(int));
40  int IfUser = 1;
41  int NElMin = 0;
42  int NElMax = 0;
43  for(int i=1;i<argc;i++){
44  if(!strcmp(*(argv+i),"--col")){
45  if(i<argc-1){
46  i+=1;
47  }
48  }
49  else if(!strcmp(*(argv+i),"--exec")){
50  if(i<argc-1){
51  sprintf(Executex,"%s",*(argv+i+1));
52  sprintf(Executey,"%s",*(argv+i+2));
53  i+=2;
54  sprintf(Comando,"exec");
55  IfUser = 0;
56  }
57  }
58  else if(!strncmp(*(argv+i),"--Border",8)){
59  if(argc < i){printf("How many values?\n");return 0;}
60  sscanf(argv[i+1],"%lf",Border);
61  sscanf(argv[i+2],"%lf",Border+1);
62  IfBorder=1;
63  i+=2;
64  }
65  else if(!strncmp(*(argv+i),"-l",2)){
66  if(argc < i){printf("Whic log scale?\n");return 0;}
67  char LogScale[12];
68  sscanf(argv[i+1],"%s",LogScale);
69  LogLog = 0;
70  if(!strcmp(LogScale,"logx"))
71  DIS_ADD_TYPE(LogLog,DIS_LOGX);
72  else if(!strcmp(LogScale,"logy"))
73  DIS_ADD_TYPE(LogLog,DIS_LOGY);
74  else if(!strcmp(LogScale,"loglog"))
75  DIS_ADD_TYPE(LogLog,DIS_LOGLOG);
76  else{
77  printf("Log scale not recognized\n");
78  }
79  i+=1;
80  }
81  else if(!strncmp(*(argv+i),"--NElMin",8)){
82  if(argc < i){printf("How many values?\n");return 0;}
83  sscanf(argv[i+1],"%d",&NElMin);
84  i+=1;
85  }
86  else if(!strncmp(*(argv+i),"--NElMax",8)){
87  if(argc < i){printf("How many values?\n");return 0;}
88  sscanf(argv[i+1],"%d",&NElMax);
89  i+=1;
90  }
91  else if(!strncmp(*(argv+i),"-m",2)){
92  if(argc < i){printf("which interval for the running average?\n");return 0;}
93  sscanf(argv[i+1],"%d",&NMedia);
94  i+=1;
95  }
96  else if(!strncmp(*(argv+i),"--norm",6)){
97  IfNorm = 1;
98  }
99  else if(!strncmp(*(argv+i),"--NSample",8)){
100  if(argc < i){printf("How many samples?\n");return 0;}
101  sscanf(argv[i+1],"%d",&NSample);
102  i+=1;
103  }
104  else if(!strncmp(*(argv+i),"-o",2)){
105  if(argc < i){printf("which output file?\n");return 0;}
106  sscanf(argv[i+1],"%s",Uscita);
107  i+=1;
108  }
109  else if(!strncmp(*(argv+i),"-p",2)){
110  if(argc < i){printf("what parameter?\n");return 0;}
111  sscanf(argv[i+1],"%lf",&Param);
112  i+=1;
113  }
114  else if(!strncmp(*(argv+i),"-v",2)){
115  if(argc < i){printf("How many values?\n");return 0;}
116  sscanf(argv[i+1],"%d",&NBin);
117  i+=1;
118  }
119  else if(!strncmp(*(argv+i),"-x",2)){
120  i+=1;
121  }
122  else if(!strncmp(*(argv+i),"-y",2)){
123  i+=1;
124  }
125  else if(!strncmp(*(argv+i),"--",2)){
126  sprintf(Comando,"%s",argv[i]+2);
127  IfUser = 0;
128  }
129  else if(strncmp(*(argv+i),"-",1)){
130  ArgFile[NFile] = i;
131  NFile++;
132  }
133  }
134  if(NFile == 0){
135  printf("File arguments are missing\n");
136  //return 0;
137  }
138  v1 = new VarDatFile(argv,ArgFile,NFile,NBin);
139  for(int i=1;i<argc;i++){
140  if(!strncmp(*(argv+i),"-a",2)){
141  v1->ImpSequence();
142  }
143  if(!strncmp(*(argv+i),"--col",5)){
144  if(argc < i){printf("Which column?\n");return 0;}
145  sscanf(argv[i+1],"%d",&CoordY);
146  v1->Punta(CoordY);
147  i+=1;
148  }
149  else if(!strncmp(*(argv+i),"-l",2)){
150  if(argc < i){printf("Which coord for the log scale?\n");return 0;}
151  char LogScale[12];
152  sscanf(argv[i+1],"%s",LogScale);
153  LogLog = 0;
154  if(!strcmp(LogScale,"logx"))
155  DIS_ADD_TYPE(LogLog,DIS_LOGX);
156  else if(!strcmp(LogScale,"logy"))
157  DIS_ADD_TYPE(LogLog,DIS_LOGY);
158  else if(!strcmp(LogScale,"loglog"))
159  DIS_ADD_TYPE(LogLog,DIS_LOGLOG);
160  else{
161  printf("Log scale not recognized\n");
162  }
163  i+=1;
164  }
165  else if(!strncmp(*(argv+i),"-y1",3)){
166  if(argc < i){printf("Which column?\n");return 0;}
167  sscanf(argv[i+1],"%d",&CoordY1);
168  i+=1;
169  }
170  else if(!strncmp(*(argv+i),"-x1",3)){
171  if(argc < i){printf("Which column?\n");return 0;}
172  sscanf(argv[i+1],"%d",&CoordX1);
173  i+=1;
174  }
175  else if(!strncmp(*(argv+i),"-x",2)){
176  if(argc < i){printf("which x coord?\n");return 0;}
177  sscanf(argv[i+1],"%d",&CoordX);
178  v1->ImpCoordX(CoordX);
179  i+=1;
180  }
181  else if(!strncmp(*(argv+i),"-y",2)){
182  if(argc < i){printf("Which column?\n");return 0;}
183  sscanf(argv[i+1],"%d",&CoordY);
184  v1->ImpCoordY(CoordY);
185  i+=1;
186  }
187  }
188  NMax = v1->pNMax();
189  if(NElMax == 0 || NElMax > NMax) NElMax = NMax;
190  if(NElMin > NMax || NElMin < 0) NElMin = 0;
191  NVar = v1->pNVar();
192  sp = new double;
193  sw = new double[NMax];
194  Intervalli = new double[NBin];
195  dInter = new double[NBin];
196  printf("Loaded a file of %d points %d cols from %d files: %s ... to write on %s\n",NMax,NVar,NFile,argv[ArgFile[0]],Uscita);
197  //goto vet;
198  v1->pMinMaxGlob(NElMin,NElMax);
199  while(strcmp(Comando,"q") ){
200  if(IfUser){
201  printf("Matematica> ");
202  scanf("%s",Comando);
203  }
204  else if(!strcmp(Comando,"agg")){
205  char NomeFile[120];
206  printf("Nome file :");
207  scanf("%s",NomeFile);
208  v1->Aggiungi(NomeFile);
209  if(!IfUser) return 0;
210  }
211  else if(!strcmp(Comando,"create")){
212  FILE *Ciccia = fopen(Uscita,"w");
213  if(Ciccia == NULL){printf("%s could not be created\n",Uscita);exit(0);}
214  int NPassi = 100;
215  double InvNPassi = 1./(double)NPassi;
216  double Freq1 = 32.;
217  double Freq2 = 4.;
218  for(int px=0;px<NPassi;px++){
219  double x = px*InvNPassi*DUE_PI;
220  fprintf(Ciccia,"%lf %lf\n",x,sin(x*Freq1));
221  }
222  }
223  else if(!strcmp(Comando,"create2d")){
224  FILE *Ciccia = fopen(Uscita,"w");
225  int NPassi = 100;
226  double InvNPassi = 1./(double)NPassi;
227  double Freq1 = 32.;
228  double Freq2 = 4.;
229  int IfOctave = 1;
230  if(!IfOctave){
231  for(int px=0;px<NPassi;px++){
232  for(int py=0;py<NPassi;py++){
233  double x = px*InvNPassi*DUE_PI;
234  double y = py*InvNPassi*DUE_PI;
235  fprintf(Ciccia,"%lf %lf %lf\n",x,y,sin(x*y*Freq1));
236  }
237  }
238  }
239  else{
240  fprintf(Ciccia,"# Created by ElPoly\n");
241  fprintf(Ciccia,"# name: tz\n");
242  fprintf(Ciccia,"# type: matrix\n");
243  fprintf(Ciccia,"# rows: %d\n",NPassi);
244  fprintf(Ciccia,"# columns: %d\n",NPassi);
245  for(int px=0;px<NPassi;px++){
246  for(int py=0;py<NPassi;py++){
247  double x = px*InvNPassi*DUE_PI;
248  double y = py*InvNPassi*DUE_PI;
249  fprintf(Ciccia,"%lf ",sin(x*y*Freq2));
250  }
251  fprintf(Ciccia,"\n");
252  }
253  }
254  fclose(Ciccia);
255  if(!IfUser) return 0;
256  }
257  else if(!strcmp(Comando,"dis")){
258  Mom = v1->DistrGaussSegnale(NElMin,NElMax,IfNorm);
259  FILE *USCITA;
260  if((USCITA = fopen(Uscita,"w"))==0){
261  printf("Non s'apre %s\n Permessi?\n",Uscita);
262  }
263  fprintf(USCITA,"# %lf %lf %lf\n",Mom.Uno,Mom.Due,Mom.Tre);
264  for(int i=0;i<NBin;i++){
265  double x = Mom.Delta*i+Mom.Min;
266  fprintf(USCITA,"%lf %g\n",x,v1->pInter(i));
267  }
268  fclose(USCITA);
269  printf("Media %lf Scarto %lf N %d\n",Mom.Uno,Mom.Due,Mom.Num);
270  if(!IfUser) return 0;
271  }
272  else if(!strcmp(Comando,"disErr")){
273  Mom = v1->DistrSignErr(NElMin,NElMax,Border,IfNorm);
274  FILE *USCITA;
275  if((USCITA = fopen(Uscita,"w"))==0){
276  printf("Non s'apre %s\n Permessi?\n",Uscita);
277  }
278  fprintf(USCITA,"# %lf %lf %lf\n",Mom.Uno,Mom.Due,Mom.Tre);
279  for(int i=0;i<NBin;i++){
280  double x = Mom.Delta*i+Mom.Min;
281  fprintf(USCITA,"%lf %g %g\n",x,v1->pInter(i),v1->pError(i));
282  }
283  fclose(USCITA);
284  printf("Media %lf Scarto %lf N %d\n",Mom.Uno,Mom.Due,Mom.Num);
285  if(!IfUser) return 0;
286  }
287  else if(!strcmp(Comando,"disB")){
288  Mom = v1->DistrSegnale(NElMin,NElMax,Border,IfNorm);
289  FILE *USCITA;
290  if((USCITA = fopen(Uscita,"w"))==0){
291  printf("Non s'apre %s\n Permessi?\n","Distribuzione.dat");
292  }
293  fprintf(USCITA,"# %lf %lf %lf\n",Mom.Uno,Mom.Due,Mom.Tre);
294  for(int i=0;i<NBin;i++){
295  double x = Mom.Delta*i+Mom.Min;
296  fprintf(USCITA,"%lf %g\n",x,v1->pInter(i));
297  }
298  fclose(USCITA);
299  printf("Media %lf Scarto %lf N %d\n",Mom.Uno,Mom.Due,Mom.Num);
300  if(!IfUser) return 0;
301  }
302  else if(!strcmp(Comando,"disLog")){
303  Mom = v1->DistrLogSegnale(NElMin,NElMax,IfNorm);
304  FILE *USCITA;
305  if((USCITA = fopen(Uscita,"w"))==0){
306  printf("Non s'apre %s\n Permessi?\n","Distribuzione.dat");
307  }
308  fprintf(USCITA,"# %lf %lf %lf\n",Mom.Uno,Mom.Due,Mom.Tre);
309  for(int i=0;i<NBin;i++){
310  double x = Mom.Delta*i+Mom.Min;
311  fprintf(USCITA,"%lf %g\n",x,v1->pInter(i));
312  }
313  fclose(USCITA);
314  printf("Media %lf Scarto %lf N %d\n",Mom.Uno,Mom.Due,Mom.Num);
315  if(!IfUser) return 0;
316  }
317  else if(!strcmp(Comando,"disExp")){
318  Mom = v1->DistrExpSegnale(NElMin,NElMax,IfNorm);
319  FILE *USCITA;
320  if((USCITA = fopen(Uscita,"w"))==0){
321  printf("Non s'apre %s\n Permessi?\n","Distribuzione.dat");
322  }
323  fprintf(USCITA,"# %lf %lf %lf\n",Mom.Uno,Mom.Due,Mom.Tre);
324  for(int i=0;i<NBin;i++){
325  double x = Mom.Delta*i+Mom.Min;
326  fprintf(USCITA,"%lf %g\n",x,v1->pInter(i));
327  }
328  fclose(USCITA);
329  printf("Media %lf Scarto %lf N %d\n",Mom.Uno,Mom.Due,Mom.Num);
330  if(!IfUser) return 0;
331  }
332  else if(!strcmp(Comando,"disSample")){
333  double xBound[2] = {Border[0],Border[1]};
334  double **Distr;
335  Distr = new double*[NSample];
336  for(int s=0;s<NSample;s++){
337  Distr[s] = new double[NBin];
338  }
339  v1->DistrSignSample(NElMin,NElMax,Distr,NSample,IfNorm,xBound);
340  FILE *USCITA;
341  if((USCITA = fopen(Uscita,"w"))==0){
342  printf("Non s'apre %s\n Permessi?\n","Distribuzione.dat");
343  }
344  fprintf(USCITA,"# distribution per sample\n");
345  for(int b=0;b<NBin;b++){
346  double x = b/(double)NBin*(xBound[1]-xBound[0]) + xBound[0];
347  fprintf(USCITA,"%lf ",x);
348  for(int s=0;s<NSample;s++){
349  fprintf(USCITA,"%lf ",Distr[s][b]);
350  }
351  fprintf(USCITA,"\n");
352  }
353  fclose(USCITA);
354  if(!IfUser) return 0;
355  }
356  else if(!strcmp(Comando,"DrawEq")){
357  FILE *FWrite = fopen(Uscita,"w");
358  if(!IfBorder){
359  Border[0] = v1->pxMinGlob(NElMin,NElMax);
360  Border[1] = v1->pxMaxGlob(NElMin,NElMax);
361  }
362  double Dx = (Border[1]-Border[0])/(double)NBin;
363  for(int i=0;i<NBin;i++){
364  double x = i*Dx+Border[0];
365  fprintf(FWrite,"%lf %lf\n",x,v1->fProva(x));
366  }
367  fclose(FWrite);
368  if(!IfUser) return 0;
369  }
370  else if(!strcmp(Comando,"coord")){
371  printf("Number of coordinate: ");
372  int NCoord = 0;
373  scanf("%d",&NCoord);
374  if(NCoord >= 0 && NCoord<NVar)
375  v1->Punta(NCoord);
376  if(!IfUser) return 0;
377  }
378  else if(!strcmp(Comando,"MatchCoord")){
379  double *st = (double *)calloc(NMax,sizeof(double));
380  double tInit = v1->Val(CoordX1,0);
381  int mOld = 0;
382  for(int n=NElMin,m=0;n<NElMax;n++){
383  mOld = m-1;
384  do {
385  double t = tInit + v1->Val(CoordX,m);
386  //printf("%d %d %lf <= %lf < %lf -> %lf\n",n,m,v1->Val(CoordX1,n),t,v1->Val(CoordX1,n+1),st[n]);
387  if(v1->Val(CoordX1,n) <= t && t < v1->Val(CoordX1,n+1))
388  break;
389  m++;
390  if(m >= NMax){
391  m = mOld;
392  break;
393  }
394  } while(1==1);
395  st[n] = v1->Val(CoordY,m);
396  //printf("%d %d %lf <= %lf < %lf -> %lf\n",n,m,v1->Val(CoordX1,n),v1->Val(CoordX,m)+tInit,v1->Val(CoordX1,n+1),st[n]);
397  m++;
398  }
399  FILE *USCITA;
400  if((USCITA = fopen(Uscita,"w"))==0)
401  printf("Non s'apre %s\n Permessi?\n",Uscita);
402  for(int i=NElMin;i<NElMax;i++){
403  fprintf(USCITA,"%lf %.5g %.5g\n",st[i],v1->Val(CoordY1,i),v1->Val(CoordY1+1,i));
404  }
405  fclose(USCITA);
406  free(st);
407  if(!IfUser) return 0;
408  }
409  else if(!strcmp(Comando,"norm")){
410  FILE *USCITA;
411  if((USCITA = fopen(Uscita,"w"))==0)
412  printf("Non s'apre %s\n Permessi?\n",Uscita);
413  v1->NormalizzaSegnale(NElMin,NElMax);
414  for(int i=NElMin;i<NElMax;i++){
415  fprintf(USCITA,"%lf %.5g\n",v1->Abscissa(CoordY,i),v1->Val(CoordY,i));
416  }
417  fclose(USCITA);
418  if(!IfUser) return 0;
419  }
420  else if(!strcmp(Comando,"sin")){
421  FILE *USCITA;
422  if((USCITA = fopen(Uscita,"w"))==0){
423  printf("Non s'apre %s\n Permessi?\n",Uscita);
424  }
425  for(int i=0;i<100;i++){
426  double x = i*DUE_PI*.01;
427  fprintf(USCITA,"%g\n",sin(x));
428  }
429  fclose(USCITA);
430  if(!IfUser) return 0;
431  }
432  else if(!strcmp(Comando,"sum")){
433  double Sum = 0.;
434  for(int i=NElMin;i<NElMax;i++){
435  Sum += v1->Val(CoordY,i);
436  }
437  printf("Sum: %lf\n",Sum);
438  if(!IfUser) return 0;
439  }
440  else if(!strcmp(Comando,"SpD")){
441  printf("Calcolo lo spettro \n");
442  m1->Spettro(sp,sw,NMax);
443  FILE *USCITA;
444  if((USCITA = fopen(Uscita,"w"))==0){
445  printf("Non s'apre %s\n Permessi?\n",Uscita);
446  return 0;
447  }
448  for(int i=0;i<NMax/2;i++){
449  fprintf(USCITA,"%g\n",sw[i]);
450  }
451  fclose(USCITA);
452  if(!IfUser) return 0;
453  }
454  else if(!strcmp(Comando,"sort")){
455  v1->Punta(0);
456  v1->Sort();
457  v1->ScriviTutto(Uscita,0,NElMin,NElMax);
458  if(!IfUser) return 0;
459  }
460  else if(!strcmp(Comando,"Cas")){
461  printf("Proviamo il generatore di numeri casuali \n");
462  FILE *USCITA;
463  if((USCITA = fopen(Uscita,"w"))==0){
464  printf("Non s'apre %s\n Permessi?\n",Uscita);
465  }
466  for(int i=0;i<10000;i++){
467  fprintf(USCITA,"%lf\n",m1->Casuale());
468  }
469  fclose(USCITA);
470  if(!IfUser) return 0;
471  }
472  else if(!strcmp(Comando,"Func")){
473  FILE *USCITA;
474  printf("Proviamo una funzione \n");
475  if((USCITA = fopen(Uscita,"w"))==0){
476  printf("Non s'apre %s\n Permessi?\n",Uscita);
477  }
478  for(double x=0.,Dx=.001,Max=5.;x<Max;x+=Dx){
479  printf("%lf\n",x/Max);
480  fprintf(USCITA,"%lf %g %g %g\n",x,m1->QuasiBessel(x,0),m1->QuasiBessel(x,1),m1->Bessel(x,0));
481  //fprintf(USCITA,"%lf %g %g %g\n",x,m1->Neumann(x,0),m1->Neumann(x,1),m1->Neumann(x,-1));
482  x += Dx;
483  }
484  fclose(USCITA);
485  if(!IfUser) return 0;
486  }
487  else if(!strcmp(Comando,"SpeLine")){
488  int NBin = NElMax-NElMin;
489  double InvNVar = .5/(double)(NVar);
490  double InvNBin = 1./(double)(NMax);
491  double *Spe = (double *)calloc(NBin,sizeof(double));
492  v1->SpeLine(NElMin,NElMax,NBin,Spe);
493  FILE *FWrite = fopen("StalkLineSpectrum.dat","w");
494  fprintf(FWrite,"# q=2\\pi n/N <h(q)^2>\n");
495  char FName[120];
496  for(int n=0;n<NBin;n++){
497  fprintf(FWrite,"%lf %g \n",DUE_PI*InvNBin*n,Spe[n]*InvNVar*InvNBin);
498  }
499  fclose(FWrite);
500  free(Spe);
501  if(!IfUser) return 0;
502  }
503  else if(!strcmp(Comando,"spe") ){
504  FILE *USCITA;
505  m1->Spettro(v1->sp,sw,NMax);
506  if((USCITA = fopen(Uscita,"w"))==0){
507  printf("Non s'apre %s\n Permessi?\n",Uscita);
508  }
509  for(int i=0;i<NMax;i++){
510  fprintf(USCITA,"%g\n",sw[i]);
511  }
512  fclose(USCITA);
513  break;
514  if(!IfUser) return 0;
515  }
516  else if(!strcmp(Comando,"RadTens") ){
517  if(NFile != 4){
518  printf("I need the three components of the radial pressure\n");
519  if(!IfUser) return 0;
520  }
521  FILE *RadTens = fopen(Uscita,"w");
522  double Vol1 = 0.;
523  double Vol2 = 0.;
524  for(int n=0;n<NMax-1;n++){
525  Vol1 = DUE_PI*SQR(v1->Abscissa(CoordY,n+1));
526  double Prr = sqrt( SQR(v1->Val(1,n)) + SQR(v1->Val(6,n)));
527  double Pzz = v1->Val(11,n);
528  double PTrace = v1->Val(1,n)+v1->Val(5,n)+v1->Val(11,n);
529  double Ptt = .5*(PTrace - Prr);
530  double SurfTens = (Prr - .5*(Ptt+Pzz))/(Vol1-Vol2);
531  fprintf(RadTens,"%lf %lf\n",v1->Abscissa(CoordY,n),SurfTens);
532  Vol2 = Vol1;
533  }
534  fclose(RadTens);
535  if(!IfUser) return 0;
536  }
537  else if(!strcmp(Comando,"perm") ){
538  int NMax = 16;
539  PERMUTE *Perm = (PERMUTE *)calloc(NMax,sizeof(PERMUTE));
540  m1->PermuteRandomAll(Perm,NMax);
541  if(!IfUser) return 0;
542  }
543  else if(!strcmp(Comando,"solve") ){
544  Matrice *A = new Matrice(4,4);
545 // Matrice Pasquale(4);
546 // Matrice Romualdo(4);
547 // Matrice Adalgiso(4);
548 // Matrice Goffredo(4);
549 // Pasquale.RandomFill(1.);
550 // Romualdo.RandomFill(1.);
551 // //Pasquale = Romualdo;
552 // Adalgiso.Mult(Pasquale,Romualdo);
553 // Adalgiso.Print();
554 // exit(0);
555  A->Set(0,0,1.);A->Set(0,1,6.);A->Set(0,2,5.);A->Set(0,3,8.);
556  A->Set(1,0,5.);A->Set(1,1,2.);A->Set(1,2,5.);A->Set(1,3,5.);
557  A->Set(2,0,6.);A->Set(2,1,3.);A->Set(2,2,7.);A->Set(2,3,1.);
558  A->Set(3,0,8.);A->Set(3,1,7.);A->Set(3,2,2.);A->Set(3,3,6.);
559 
560 // A->Set(0,0,1.);A->Set(0,1,0.);A->Set(0,2,0.);A->Set(0,3,1.);
561 // A->Set(1,0,0.);A->Set(1,1,1.);A->Set(1,2,0.);A->Set(1,3,0.);
562 // A->Set(2,0,0.);A->Set(2,1,0.);A->Set(2,2,1.);A->Set(2,3,0.);
563 // A->Set(3,0,2.);A->Set(3,1,0.);A->Set(3,2,0.);A->Set(3,3,1.);
564  // A->RandomFill(20.);
565  A->Print();
566  double *Known = (double *)calloc(4,sizeof(double));
567  Known[0] = 1.; Known[1] = 4.; Known[2] = 3.; Known[3] = 2.;
568  double *UnKnown = (double *)calloc(4,sizeof(double));
569  A->Solve(Known,UnKnown);
570  for(int i=0;i<4;i++){
571  printf("%lf %lf\n",Known[i],UnKnown[i]);
572  }
573  if(!IfUser) return 0;
574  }
575  else if(!strcmp(Comando,"zeri") ){
576  double a=0.,b=0.;
577  int rad;
578 // printf("Minimo valore: ");
579 // scanf("%lf",&a);
580 // printf("Massimo valore: ");
581 // scanf("%lf",&b);
582 //%% 2.5
583 //%% 2.0
584 //%% 1.5 23.376860 3817.700000
585 //%% 1.0 21.643605 3817.700000
586  int NRadici = 4;
587  double *Radici;
588  m1->Ypsilon = 23.37;
589  m1->PreFact = 3817.7;;
591  Radici = (double *)malloc(NRadici*sizeof(double));
592  rad = m1->Zeri(0.,DUE_PI/2.,Radici,NRadici);
593  if(!IfUser) return 0;
594  }
595  else if(!strcmp(Comando,"fal") ){
596  double a=0.,b=0.;
597  RADICE Zero;
598 // printf("Minimo valore: ");
599 // scanf("%lf",&a);
600 // printf("Massimo valore: ");
601 // scanf("%lf",&b);
602  Zero = m1->RegulaFalsi(.1,3.);
603  printf("Lo Zero e in x %lf\n",Zero.Zero);
604  if(!IfUser) return 0;
605  }
606  else if(!strcmp(Comando,"new") ){
607  double a=0.;
608  RADICE Zero;
609  printf("Punto iniziale: ");
610  scanf("%lf",&a);
611  Zero = m1->Newton(a);
612  printf("Lo Zero e in x %lf\n",Zero.Zero);
613  if(!IfUser) return 0;
614  }
615  else if(!strcmp(Comando,"fun") ){
616  FILE *USCITA;
617  USCITA = fopen("Funzione.dat","w");
618  double A = 0.;
619  double B = DUE_PI;
620  int NPunti = 10000;
621  double Delta = (B-A)/(double)NPunti;
622  m1->Ypsilon = 3.814;
623  m1->PreFact = 0.169370;
625  for(int i=0;i<NPunti;i++){
626  fprintf(USCITA,"%g %g\n",Delta*(double)i,m1->Evalx(Delta*(double)i));
627  }
628  fclose(USCITA);
629  if(!IfUser) return 0;
630  }
631  else if(!strcmp(Comando,"FitLin") ){
632  FILE *USCITA = fopen(Uscita,"w");
633  RETTA r1 = v1->InterRettSegnale(CoordY,NElMin,NElMax,LogLog);
634  printf("y = (%lf +- %lf)x + (%lf +- %lf)\n",r1.m,r1.ErrM,r1.q,r1.ErrQ);
635  fprintf(USCITA,"#y = (%lf +- %lf)x + (%lf +- %lf)\n",r1.m,r1.ErrM,r1.q,r1.ErrQ);
636  if(!LogLog){
637  for(int i=0;i<NMax;i++){
638  fprintf(USCITA,"%lf %lf %lf\n",v1->Abscissa(CoordY,i),v1->Val(CoordY,i),v1->Abscissa(CoordY,i)*r1.m+r1.q);
639  }
640  }
641  else{
642  for(int i=0;i<NMax;i++){
643  if(v1->Abscissa(CoordY,i) <= 0.) continue;
644  if(v1->Val(CoordY,i) <= 0.) continue;
645  fprintf(USCITA,"%lf %lf %lf\n",log10(v1->Abscissa(CoordY,i)),log10(v1->Val(CoordY,i)),log10(v1->Abscissa(CoordY,i))*r1.m+r1.q);
646  }
647  }
648  fclose(USCITA);
649  if(!IfUser) return 0;
650  }
651  else if(!strcmp(Comando,"FitPara") ){
652  FILE *USCITA = fopen(Uscita,"w");
653  PARABOLA p1 = v1->ParabolaSegnale(CoordY,NElMin,NElMax,LogLog);
654  printf("%.2g + %.2g x + %.2g x^2\n",p1.a0,p1.a1,p1.a2);
655  printf("Minimum %.2g \n",p1.Minimo);
656  fprintf(USCITA,"#%.2g + %.2g x + %.2g x^2\n",p1.a0,p1.a1,p1.a2);
657  fprintf(USCITA,"#Minimum %.2g \n",p1.Minimo);
658  for(int i=0;i<NMax;i++){
659  fprintf(USCITA,"%lf %lf %lf\n",v1->Abscissa(CoordY,i),v1->Val(CoordY,i),v1->Abscissa(CoordY,i)*p1.a1+v1->Abscissa(CoordY,i)*v1->Abscissa(CoordY,i)*p1.a2+p1.a0);
660  }
661  fclose(USCITA);
662  if(!IfUser) return 0;
663  }
664  else if(!strcmp(Comando,"FitExp") ){
665  FILE *USCITA = fopen(Uscita,"w");
666  RETTA r1 = v1->InterExpSegnale(CoordY,NElMin,NElMax,LogLog);
667  printf("y=(%.3g+-%.3g)exp(x/(%.3g+-%.3g))\n",log10(r1.q),log10(r1.ErrQ),1./r1.m,r1.ErrM/SQR(r1.m));
668  printf("decay t_d = %.4g +- %.4g\n",1./r1.m,r1.ErrM/SQR(r1.m));
669  fprintf(USCITA,"#y=(%.3g+-%.3g)exp(x/(%.3g+-%.3g))\n",log10(r1.q),log10(r1.ErrQ),1./r1.m,r1.ErrM/SQR(r1.m));
670  fprintf(USCITA,"#decay t_d = %.4g +- %.4g\n",1./r1.m,r1.ErrM/SQR(r1.m));
671  for(int i=0;i<NMax;i++){
672  fprintf(USCITA,"%lf %lf %lf\n",v1->Abscissa(CoordY,i),v1->Val(CoordY,i),exp(r1.q+r1.m*v1->Abscissa(CoordY,i)));
673  }
674  fclose(USCITA);
675  if(!IfUser) return 0;
676  }
677  else if(!strcmp(Comando,"Prob2Nrg") ){
678  FILE *File2Write = fopen(Uscita,"w");
679  if(File2Write == NULL){
680  printf("Could not open %s\n",Uscita);
681  return 1;
682  }
683  MOMENTI Mom = v1->DistrSegnale(NElMin,NElMax,IfNorm);
684  v1->NormalizzaInter();
685  for(int v=0;v<NBin;v++){
686  double x = Mom.Min + v*Mom.Delta;
687  if(v1->pInter(v) > 0.)
688  fprintf(File2Write,"%lf %lf\n",x,-log(v1->pInter(v)));
689  }
690  fclose(File2Write);
691  if(!IfUser) return 0;
692  }
693  else if(!strcmp(Comando,"min") ){
694  double Min = v1->Val(CoordY,0);
695  int nMin = 0;
696  for(int n=0;n<NMax;n++){
697  if(Min > v1->Val(CoordY,n)){
698  Min = v1->Val(CoordY,n);
699  nMin = n;
700  }
701  }
702  printf("%lf %lf\n",v1->Val(CoordX,nMin),v1->Val(CoordY,nMin));
703  if(!IfUser) return 0;
704  }
705  else if(!strcmp(Comando,"mm") ){
706  printf("Number of points per sample: %d\n",NMedia);
707  int NSample = v1->MediaMobSegnale(NMedia);
708  double Deltax = 1.;//(v1->pxMax()-v1->pxMin())/(double)NSample;
709  FILE *USCITA;
710  if((USCITA = fopen(Uscita,"w"))==0)
711  printf("Non s'apre %s\n Permessi?\n",Uscita);
712  for(int i=0;i<NSample;i++){
713  fprintf(USCITA,"%lf %lf %lf\n",v1->Abscissa(CoordY,i*NMedia)+Deltax,v1->pPunti(i),v1->pPuntiErr(i));
714  }
715  fclose(USCITA);
716  if(!IfUser) return 0;
717  }
718  else if(!strcmp(Comando,"WeightHisto") ){
719  v1->WeightHistoSign(NFile);
720  // FILE *USCITA;
721  // if((USCITA = fopen(Uscita,"w"))==0)
722  // printf("Non s'apre %s\n Permessi?\n",Uscita);
723  // for(int i=0;i<NSample;i++){
724  // fprintf(USCITA,"%lf %lf %lf\n",v1->Abscissa(CoordY,i*NMedia)+Deltax,v1->pPunti(i),v1->pPuntiErr(i));
725  // }
726  // fclose(USCITA);
727  if(!IfUser) return 0;
728  }
729  else if(!strcmp(Comando,"int")){
730  double Int = v1->IntSegnale();
731  FILE *USCITA;
732  if((USCITA = fopen(Uscita,"w"))==0)
733  printf("Non s'apre %s\n Permessi?\n",Uscita);
734  printf("%d %d %d\n",NElMin,NElMax,v1->pNMax());
735  for(int i=NElMin;i<NElMax;i++){
736  fprintf(USCITA,"%lf %lf\n",v1->Abscissa(CoordY,i),v1->pPunti(i));
737  }
738  printf("%lf\n",Int);
739  if(!IfUser) return 0;
740  }
741  else if(!strcmp(Comando,"mom")){
742  for(int v=0;v<NVar;v++){
743  v1->Punta(v);
744  Mom = v1->DistrGaussSegnale(NElMin,NElMax,IfNorm);
745  printf("%d) %lf pm %lf\n",v,Mom.Uno,Mom.Due);
746  }
747  if(!IfUser) return 0;
748  }
749  else if(!strcmp(Comando,"exec")){
750  if(IfUser){
751  printf("Inserire la formula per x: ");
752  scanf("%s",Executex);
753  printf("Inserire la formula per y: ");
754  scanf("%s",Executey);
755  }
756  FILE *Ciccia = fopen(Uscita,"w");
757  for(int i=NElMin;i<NElMax;i++){
758  for(int v=0;v<NVar;v++){
759  double Val = v1->Val(v,i);
760  if(v==CoordX)
761  Val = v1->ExecFormula(v1->Val(CoordX,i),v1->Val(CoordY,i),Executex);
762  else if(v==CoordY)
763  Val = v1->ExecFormula(v1->Val(CoordX,i),v1->Val(CoordY,i),Executey);
764  fprintf(Ciccia,"%.5g ",Val);
765  }
766  fprintf(Ciccia,"\n");
767  }
768  fclose(Ciccia);
769  if(!IfUser) return 0;
770  }
771  else if(!strcmp(Comando,"widom")){
772  FILE *Widom = fopen(Uscita,"w");
773  for(int i=NElMin;i<NElMax;i++){
774  if(v1->Val(3,i)*v1->Val(1,i) > 0.){
775  double Arg = v1->Val(3,i)/(v1->Val(1,i));
776  double ChemPot = v1->Val(0,i) + log(Arg);
777  fprintf(Widom,"%lf %lf\n",v1->Val(0,i),ChemPot);
778  }
779  }
780  fclose(Widom);
781  if(!IfUser) return 0;
782  }
783  else if(!strcmp(Comando,"DecayTime")){
784  double *st = new double[NVar];
785  double *sw = new double[NVar];
786  FILE *DecayTime = fopen(Uscita,"w");
787  for(int f=0;f<NVar;f++){
788  for(int i=NElMin;i<NElMax;i++){
789  st[f] = v1->Val(0,i)*v1->Val(f,i);
790  }
791  st[f] /= (double)v1->pNMax();
792  }
793  m1->Spettro(st,sw,NVar);
794  for(int f=0;f<NVar;f++){
795  fprintf(DecayTime,"%d %g %g\n",f,sw[f],st[f]);
796  }
797  fclose(DecayTime);
798  if(!IfUser) return 0;
799  }
800  else if(!strcmp(Comando,"reverse") ){
801  v1->Reverse();
802  v1->ScriviTutto(Uscita,0,NElMin,NElMax);
803  if(!IfUser) return 0;
804  }
805  else if(!strcmp(Comando,"DoubleDist") ){
806  v1->DoubleDistFluct();
807  v1->ScriviPunti(Uscita);
808  if(!IfUser) return 0;
809  }
810  else if(!strcmp(Comando,"smooth") ){
811  v1->Smooth(Param,CoordY,NElMin,NElMax);
812  v1->Smooth(Param,CoordY,NElMin,NElMax);
813  v1->Smooth(Param,CoordY,NElMin,NElMax);
814  v1->ScriviPunti(Uscita);
815  if(!IfUser) return 0;
816  }
817  else if(!strcmp(Comando,"smoothGauss") ){
818  v1->SmoothGauss(Param,CoordY,NElMin,NElMax);
819  v1->ScriviTutto(Uscita,0,NElMin,NElMax);
820  if(!IfUser) return 0;
821  }
822  else if(!strcmp(Comando,"salva") ){
823  v1->ScriviTutto(Uscita,0,NElMin,NElMax);
824  if(!IfUser) return 0;
825  }
826  else if(!strcmp(Comando,"tecplot") ){
827  v1->TecPlot(Uscita);
828  if(!IfUser) return 0;
829  }
830  else if(!strcmp(Comando,"txvl") ){
831  v1->ExportTxvl(Uscita,NElMin,NElMax);
832  if(!IfUser) return 0;
833  }
834  else if(!strcmp(Comando,"ResBulk") ){
835  v1->RescaleToBulk(Uscita);
836  if(!IfUser) return 0;
837  }
838  else if(!strcmp(Comando,"conv") ){
839  FILE *FWrite = fopen(Uscita,"w");
840  double x = 0.;
841  for(int n=0;n<NMax;n++){
842  x += .1;
843  fprintf(FWrite,"%lf %lf\n",x,v1->Val(0,n)*1e12);
844  }
845  fclose(FWrite);
846  if(!IfUser) return 0;
847  }
848  else if(!strcmp(Comando,"salvacol") ){
849  v1->ScriviFile(Uscita,CoordY,0,NElMin,NElMax);
850  if(!IfUser) return 0;
851  }
852  else if(!strcmp(Comando,"col") ){
853  printf("Number of column: ");
854  int NCol=0;
855  scanf("%d",&NCol);
856  v1->Punta(NCol);
857  if(!IfUser) return 0;
858  }
859  else if(!strcmp(Comando,"print") ){
860  v1->Print();
861  if(!IfUser) return 0;
862  }
863  else if(!strcmp(Comando,"vet")){
864  Vettore v(1.,0.,1.);
865  Vettore u(0.,1.,0.);
866  Vettore w(0.,0.,1.);
867  Vettore n(3);
868  v.VetV(&u);
869  //printf("%lf %lf %lf\n",n.Angle(v,u)*RAD_DEG,n.SinAngle(v,u),n.CosAngle(v,u));
870  n.Print();
871  break;
872  if(!IfUser) return 0;
873  }
874  else if(!strcmp(Comando,"OrdAv")){
875  double *Distr = new double[NBin];
876  if(!IfBorder){
877  v1->pMinMaxGlob(NElMin,NElMax);
878  v1->pGlobBorder(Border,Border+1,Border+2,Border+3);
879  printf("%lf %lf\n",Border[0],Border[1]);
880  Border[1] = Border[1] - Border[0];
881  Border[1] = 1./(1.001*Border[1]);
882  }
883  else{
884  Border[1] = Border[1] - Border[0];
885  Border[1] = 1./(1.001*Border[1]);
886  }
887  printf("%lf %lf\n",Border[0],Border[1]);
888  v1->AverageOrdinate(NElMin,NElMax,Distr,Border);
889  FILE *WFile = fopen(Uscita,"w");
890  for(int b=0;b<NBin;b++){
891  double x = b/(Border[1]*NBin) + Border[0];
892  fprintf(WFile,"%lf %lf\n",x,Distr[b]);
893  }
894  free(Distr);
895  fclose(WFile);
896  if(!IfUser) return 0;
897  }
898  else if(!strcmp(Comando,"WAv")){
899  MOMENTI m1 = v1->WeightAverageSet(CoordY,NElMin,NElMax);
900  printf("Mean %lf +- %lf\n",m1.Uno,m1.Due);
901  if(!IfUser) return 0;
902  }
903  else if(!strcmp(Comando,"q") ){}
904  else{
905  printf("Comando sbagliato, riprovare\n");
906  IfUser = 1;
907  }
908  }
909  printf("Te se qe ve be te ne?\n");
910  delete m1;
911  return 0;
912 }
void Sort()
Sort with respect to a coordinate.
Definition: VarDatFile.cpp:733
double pPunti(int n)
Value of the elaborated array at position.
Definition: VarDatFile.cpp:673
double a1
a0 + a1 x + a2 x^2
void ScriviPunti(char *file)
Writes the content of the elaborated array.
Definition: VarDatFile.cpp:745
double Ypsilon
External parameter to calculate the contact angle.
Definition: Matematica.h:330
int pNMax()
Print the maximum number of data in the array.
Definition: VarDatFile.h:205
int Solve(double *Known, double *UnKnown)
Solve a system A|b = y.
void AverageOrdinate(int ElMin, int ElMax, double *Distr, double *xBound)
Average of all ordinate.
Definition: VarDatFile.cpp:325
Geometrical operations on vectors.
Definition: MatematicaVect.h:9
int Zeri(double a, double b, double *Radici, int NRadici)
Find the.
double Min
Minimum value.
void Print()
Prints the content in memory.
Definition: VarDatFile.cpp:719
Parabolas coefficients.
double Due
Second moment.
double pPuntiErr(int n)
Value of the error array at postion.
Definition: VarDatFile.cpp:680
double Tre
Third moment.
STL namespace.
double a2
a0 + a1 x + a2 x^2
RETTA InterRettSegnale(int CoordY, int NElMin, int NElMax, int LogLog)
Computes a liner fit.
Definition: VarDatFile.cpp:500
Indices for the permutation.
double fProva(double x)
Trial function.
void Punta(int n)
Points to a different column.
Definition: VarDatFile.cpp:638
void SpeLine(int ElMin, int ElMax, int NBin, double *Spe)
Spectrum of a sequence of lines.
Definition: VarDatFile.cpp:394
int PermuteRandomAll(int *Sequence, int NMass)
Permutes a sequence without repeating.
int WeightHistoSign(int NHisto)
Weighted histogram analysis.
Definition: VarDatFile.cpp:578
double * sp
Points to a column of st.
Definition: VarDatFile.h:291
double m
y = m*x + q
double QuasiBessel(double Val, int Ord)
A faster Bessel.
double Bessel(double Val, int Ord)
Bessel function.
double ErrQ
Error on the intercept.
MOMENTI WeightAverageSet(int CoordY, int ElMin, int ElMax)
Weighted distribution.
Definition: VarDatFile.cpp:371
FUNC Func
Pointer to a function.
Definition: Matematica.h:111
Linear interpolation.
Moments of a distribution.
void pGlobBorder(double *xMin, double *xMax, double *yMin, double *yMax)
Print global borders.
Definition: VarDatFile.h:209
double Zero
Point of the zero.
void ImpSequence(int Col)
Set the natural sequence as abscissa for the column Col.
Definition: VarDatFile.cpp:912
void ImpCoordY(int Ext)
Set the column for the y array.
Definition: VarDatFile.cpp:938
void DoubleDistFluct()
Distribution of distances.
double Casuale()
Random uniform number.
double ErrM
Error on the slope.
PARABOLA ParabolaSegnale(int CoordY, int NElMin, int NElMax, int LogLog)
Computes a parabolic fit.
Definition: VarDatFile.cpp:541
double ContactAngle(double x)
Definition of the contact angle.
void DistrSignSample(int ElMin, int ElMax, double **Distr, int NSample, int IfNorm, double *xBound)
Distributions of data organized in different samples.
Definition: VarDatFile.cpp:365
double IntSegnale()
Numerical integration of the signal.
Definition: VarDatFile.cpp:427
double pInter(int n)
Print the postion.
Definition: VarDatFile.h:246
double PreFact
External parameter in the definition of the contact angle.
Definition: Matematica.h:332
void ImpCoordX(int vAbs)
Set the columns for the x array.
Definition: VarDatFile.cpp:929
RADICE Newton(double a)
Use Newton to find the roots.
void RescaleToBulk(char *FName)
Rescale to bulk.
Definition: VarDatFile.cpp:788
MOMENTI DistrExpSegnale(int NElMin, int NElMax, int IfNorm)
Exponential distribution of a signal.
Definition: VarDatFile.cpp:352
MOMENTI DistrGaussSegnale(int ElMin, int ElMax, int IfNorm)
Tests a gaussian distribution of the signal.
Definition: VarDatFile.cpp:368
double a0
a0 + a1 x + a2 x^2
int pNVar()
Print the number of variables.
Definition: VarDatFile.h:203
void Reverse()
Reverse the sets.
double pError(int n)
Print the postion.
Definition: VarDatFile.h:248
double Val(int CoordY, int n)
Value at the position.
Definition: VarDatFile.cpp:657
void ScriviFile(char *file, int CoordY, int LogLog, int NVisMin, int NVisMax)
Writes the content of the pointed vector.
Definition: VarDatFile.cpp:753
Reads and stores a data file to be elaborated via Matematica.
Definition: VarDatFile.h:26
Matrice computes the algebric operations on matrices.
bool Set(int row, int column, double Val)
Set a coefficient.
Where a root was searched.
void ExportTxvl(char *FName, int NElMin, int NElMax)
Export in the txvl file format.
Definition: VarDatFile.cpp:820
double Minimo
Minimum of the parabola.
double Uno
First moment.
void SmoothGauss(double Fact, int CoordY, int NVisMin, int NVisMax)
Smooth the line.
int NormalizzaInter()
Normalize the intervals.
Definition: VarDatFile.cpp:415
double Evalx(double x)
Pointer to a generic function.
Definition: Matematica.h:113
double ExecFormula(double x, double y, char *cmd)
Execute a formula.
double pxMinGlob(int NVisMin, int NVisMax)
Print the abscissa minimum of all sets.
int Smooth(double Fact, int CoordY, int NVisMin, int NVisMax)
Smooth the line.
void TecPlot(char *FName)
Export a contour plot.
Definition: VarDatFile.cpp:859
Implementation of useful algorythms.
Definition: Matematica.h:76
double pxMaxGlob(int NVisMin, int NVisMax)
Print the abscissa maximum of all sets.
int Num
Number of points considered.
MOMENTI DistrLogSegnale(int NElMin, int NElMax, int IfNorm)
Logarithmic distribution of a signal.
Definition: VarDatFile.cpp:345
int Aggiungi(char *file)
Add the content of the file into the memory.
Definition: VarDatFile.cpp:189
void pMinMaxGlob(int NVisMin, int NVisMax)
Find the borders globally.
Definition: VarDatFile.cpp:986
RADICE RegulaFalsi(double a, double b)
Use regula falsi algorithm to find the roots.
double Abscissa(int CoordY, int n)
Value of the Abscissa at position.
Definition: VarDatFile.cpp:667
int NormalizzaSegnale(int NElMin, int NElMax)
Normalize the signal.
Definition: VarDatFile.cpp:418
RETTA InterExpSegnale(int CoordY, int NElMin, int NElMax, int LogLog)
Computes a exponential fit.
Definition: VarDatFile.cpp:530
void Print()
Print the entries.
double q
y = m*x + q
MOMENTI DistrSignErr(int NElMin, int NElMax, double *Border, int IfNorm)
Distribution of a signal with error.
Definition: VarDatFile.cpp:362
double Delta
Step interval.
void Spettro(double *st, double *sw, int NMass)
Compute the spectrum.
MOMENTI DistrSegnale(int NElMin, int NElMax, int IfNorm)
Distribution of a signal.
Definition: VarDatFile.cpp:342
void ScriviTutto(char *file, int LogLog, int NVisMin, int NVisMax)
Writes all the content in memory.
Definition: VarDatFile.cpp:766
int MediaMobSegnale(int)
Running average.
Definition: VarDatFile.cpp:573