Matlab and awk: cross section calculations etc
Matlab: prepare test data for simpson or similar
------------------------------------------------
z=ones(31,3)
z(:,2)=(0:pi/30:pi)'
z(:,3)=sin((0:pi/30:pi)')
save 'z.t' z -ascii
Matlab: Make a file with interpolated values:
------------------------------------------------
%
% x=[ 0 0.2 0.3 0.4 0.5 0.6 0.9 1.0 1.3 ];
% y=sin(2*x);
% plot (x,y,'-')
% figure
% xi=0:0.02:1.3;
% yi = interp1(x,y,xi);
% plot (xi,yi,'-')
% yi = interp1(x,y,xi,'spline');
load res;
[rows cols]=size(res);
xmax=res(rows,1);
xmin=res(1,1);
X=( xmin:(xmax-xmin)/(10*rows):xmax)';
x= res(:,1);
[nrows dum]=size(X)
nres=ones(nrows,cols);
nres(:,1)=X;
y=x;
for ui=2:cols
y=res(:,ui);
Y=interp1(x,y,X,'spline');
nres(:,ui)=Y;
end
save 'res.t' nres -ascii
------------------------------------------------
Some awk uses
------------------------------------------------
If you have the list of files:
prepare a script (see below)
------------------------------
awk ' {if (1==1) printf( "cross1 %s \n",$0); } ' THE_list_of_files
list all the files starting with d
prepare a 'script' which will use program cross1
with command line the name of the file
pipe it to csh
--------------------------------------------------
ls d* | awk ' {if (1==1) printf( "cross1 %s \n",$0); } ' | csh
Print the Last line before the line containing "velo"
-----------------------------------------------------
awk '1==1 {SEAR="velo"} $1==SEAR || $2==SEAR {print LAST} 1==1{LAST=$0}' testf
--------------------------------------------------------------
C-program cross1.c
-------------------------
#include
FILE *fp;
char filename[80];
int i,j;
float res[ 30 ],fun1[ 30 ][ 7 ],k,l,number;
main(argc,argv)
int argc; char *argv[];
{
if(argc<2)
exit();
else
{
sprintf(filename,"%s",argv[1]);
fp = fopen(filename,"r");
for (i = 0;i <= 29;i++)
{
for (j = 0;j <= 6;j++)
{
fscanf(fp,"%f",&number);
fun1[ i ][ j ] = number ;
}
}
fclose(fp); for
(i = 0;i <= 29;i++)
{
for (j = 0;j <= 6;j++)
{
printf("%f ",fun1[ i ][ j ] );
}
printf("\n");
}
}
Simpson rule integration
------------------------------------------
#include
FILE *fp, *tmpfp;
char filename[80], Command[80];
int i,j,n,m,col, OfLines, OfWords , OfChars;
float fun1[ 50 ][ 10 ],k,l,b1,b2,I1,I2,number,tcross,totcross,deltaB;
main()
{
printf("enter the file name\n");
scanf("%s",filename);
sprintf(Command,"wc %s>TMPwc",filename);
system(Command);
tmpfp=fopen("TMPwc","r");
fscanf(tmpfp,"%d %d %d", &OfLines, &OfWords , &OfChars);
n=OfLines;
if(OfWords % OfLines != 0 )
{ printf("wrong data in file %s\n", filename);
exit();
}
fclose(tmpfp);
m=OfWords / OfLines;
/*scanf("%d %d",&n,&m);*/
printf(" enter the number of the column of interest col= \n");
printf(" Note: the fist column has the number 0 \n");
scanf("%d",&col);
fp = fopen(filename,"r");
for (i = 0;i <= n-1;i++)
{
for (j = 0;j <= m-1;j++)
{
fscanf(fp,"%f",&number);
fun1[ i ][ j ] = number ;
}
}
deltaB=fun1[1][ 1]-fun1[0][ 1];
for( j=2;j0.001*deltaB) ||
((fun1[j][ 1]-fun1[j-1][ 1] - deltaB)<(-0.001*deltaB)) )
{printf("wrong imp.params in file %s\n", filename);
exit(); }
}
fclose(fp);
totcross = 0.0;
for (i = 0; i< n;i++) totcross += fun1[ i ][ col ];
for (i = 1; i< n-1;i++) totcross += fun1[ i ][ col ];
for (i = 1; i< n-1;i+=2) totcross += 2.0*fun1[ i ][ col ];
totcross *= (deltaB/3.0);
totcross = 2*3.14159*totcross;
printf("The total cross section for the event \n");
printf("described by the column %d is :\n",col);
printf("%f \n",totcross);
}
---------------------------------------------------------
In-files and out-files as command line:
---------------------------------------------------------
/*
This program accepts the name of the input file
as a command line parmeter
*/
#include
#define True (-1)
#define False (0)
FILE *fp, *outfp;
int i,j, PrintFlag;
float res[ 30 ],fun1[ 30 ][ 7 ],k,l,number;
main(argc,argv)
int argc; char *argv[];
{
if(argc<2)
exit();
else
{
if(argc==3) PrintFlag=True; else PrintFlag=False;
if(PrintFlag) outfp=fopen(argv[2],"w");
fp = fopen(argv[1],"r");
for (i = 0;i <= 29;i++)
{
for (j = 0;j <= 6;j++)
{
fscanf(fp,"%f",&number);
fun1[ i ][ j ] = number ;
}
}
fclose(fp);
for (j = 0;j <= 29;j++)
{
k = fun1[ j ][ 1 ];
l = fun1[ j ][ 2 ];
res[ j ]= k*l;
if(PrintFlag)
fprintf(outfp,"%f \n",res[ j ]);
else
printf("%f \n",res[ j ]);
}
for (j = 0;j <= 29;j++)
if(PrintFlag)
fprintf(outfp,"%f \n",fun1[j][ 3 ]);
else
printf("%f \n",fun1[j][ 3 ]);
}
}