Main Page | Alphabetical List | Data Structures | Directories | File List | Data Fields | Globals | Related Pages

libcpixe.h

Go to the documentation of this file.
00001 /***************************************************************************
00002  *   Copyright (C) 2005 by Carlos Pascual-Izarra & Miguel A. Reis          *
00003  *   carlos.pascual@itn.pt                                                 *
00004  *                                                                         *
00005  *   This program is free software; you can redistribute it and/or modify  *
00006  *   it under the terms of the GNU General Public License as published by  *
00007  *   the Free Software Foundation; either version 2 of the License, or     *
00008  *   (at your option) any later version.                                   *
00009  *                                                                         *
00010  *   This program is distributed in the hope that it will be useful,       *
00011  *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
00012  *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
00013  *   GNU General Public License for more details.                          *
00014  *                                                                         *
00015  *   You should have received a copy of the GNU General Public License     *
00016  *   along with this program; if not, write to the                         *
00017  *   Free Software Foundation, Inc.,                                       *
00018  *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
00019  ***************************************************************************/
00020 
00021 /**\file libcpixe.h
00022 libcpixe.h is the header file for the LibCPIXE library.
00023 All constants, structs, and functions are declared here.
00024 This must be included by any program using libcpixe.c
00025 */
00026 
00027 /** @todo Modify fscanf functions to check for length when reading filenames of fixed maximum length
00028 */
00029 
00030 #include "compilopt.h"
00031 
00032 #define LONGSTRINGLENGTH 256
00033 #define FILENMLENGTH LONGSTRINGLENGTH
00034 
00035 //conversion factor from degrees to rads
00036 #define DEG2RAD (M_PI/180.0)
00037 //Energy (in MeV) corresponding to 1 amu through E=mc2
00038 #define MEVAMU 931.52181 
00039 #define UATMASS (1./1822.887)
00040 
00041 //Barn to m2 conversion
00042 #define BARNTOM2  1.0e-28
00043 
00044 
00045 typedef char STATFILENAME[FILENMLENGTH] ; ///<Fixed length string for use with filenames
00046 
00047 typedef char ChemSymb[3];///<Fixed length string for use with chem. symbols
00048 
00049 typedef int atomicnumber;///<Convenient definition of atomic numbers as ints
00050 
00051 //typedef double tenarray[11];
00052 
00053 //typedef double ninearray[10];
00054 
00055 //typedef double threearray[4];
00056 
00057 //typedef double LArray[4][4];
00058 
00059 //typedef char LineNmType[6][4][11];
00060 
00061 typedef double SecXL[4];
00062 
00063 typedef enum { K, L, M } FwTipo;
00064 
00065 /*Structures*/
00066 
00067 /**Chemical element*/
00068 typedef struct ELEMENT
00069   {
00070   char symbol[3]; ///<Chemical symbol
00071   int Z;          ///<Atomic number
00072   int A;          ///<Number of nucleons
00073   double M;       ///<Natural Atomic weight, in amu 
00074   double IM;      ///<Isotopic mass, in amu
00075   }ELEMENT;
00076 
00077 /**Compound formed by various elements*/
00078 typedef struct COMPOUND
00079   {
00080   char name[30];  /**<Name of the compound*/
00081   double sumX;    /**<sum of all ammounts in the compound (useful for normalizing)*/
00082   int nelem;      /**<number of elements present in the compound*/
00083   ELEMENT *elem;  /**<Pointer to an array of elements*/
00084   double *X;      /**<Pointer to an array containing the atomic ammount of each element*/
00085   double *xn;     /**<Pointer to an array containing the NORMALIZED atomic ammount of each element*/
00086   double *w;      /**<Pointer to an array containing the NORMALIZED MASS ammount of each element**/
00087   }COMPOUND;
00088 
00089 /**Stopping Power Table.
00090   This is a modified version from that in hotstop.h*/
00091 typedef struct SPT
00092   {
00093   char ID[30];
00094   int nrows;    ///<number of rows of the table (dimension)
00095   int logmode;  ///<flag. 0 if constant Energy step (linear) and 1 if Energy is logarithmic 
00096   double Emax;  ///<Maximum energy (minimum is always 0)
00097   double Estep; ///<Step between energies, in keV (log step in case of logmode=1)
00098   double Smax;  ///<Maximum valu of the stopping power in the region E=(0,Emax)
00099   double *E;    ///<CAUTION: E will be assumed to have following units: keV
00100   double *S;    ///<CAUTION: S will be assumed to have following units: keV/1e15atcm-2
00101   double *dSdE; ///< derivative for linear interpolations: dSdE[i]=(S[i+1]-S[i])/(E[i+1]-E[i]) 
00102   }SPT;
00103 
00104 /**Filter description
00105 @todo document this structure*/
00106 typedef struct absorber
00107   {
00108   int number;
00109   char name[6];
00110   int nreg;
00111   double ereg[11];
00112   double coef[11][3];   
00113   double h[3];
00114   } absorber;
00115 
00116 
00117 /**Simulation parameters*/
00118 typedef struct SIM_PARAM
00119   {
00120   int MaxZinsample;  ///<larger Atomic number present in the whole sample
00121   int AllowSXFCorr;   ///<Flag for secondary fluorescence calculation.1=true, 0=false
00122   int AllowXEqCalc;    ///<Flag for X-ray equivalent calculation. 1=true, 0=false
00123   double DTCC;        ///<Dead Time Correction coefficient
00124   //double DTCCerr;     ///<Error in Dead Time Correction coefficient (Not really used)
00125   double ColCharge; ///<Collected charge (fluence*charge_state), in uC.
00126   double CalEner;   ///<Energy at which the calibration was done.
00127   int useFilter; ///<Flag: whether a filter on the detector is used (1) or not (0).
00128   }SIM_PARAM;
00129 
00130   
00131       
00132 /**Experimental parameters*/
00133 typedef struct EXP_PARAM
00134   {
00135   ELEMENT ion;      ///<Beam type @todo For the moment only H ions are supported.
00136   double BeamEner;  ///<Beam Energy, in keV
00137   double BeamCol;   ///<colimator aperture (diameter of the beam), in mm
00138   double DetColFac; ///<Relative diameter of det colimator to that used during calib
00139   double IncAng;    ///<Incident angle (beam to surface normal), in degrees 
00140   double DetAng;    ///<Detector Angle (detector to surface normal), in degrees
00141   double cosDet;    ///<For optimization. Cos(IncDet)
00142   double cosInc;    ///<For optimization. Cos(IncAng)
00143   double cosFac;    ///<For optimization. Cos(IncAng)/Cos(IncDet)
00144   double BeamCross; ///<Beam cross section, in cm2
00145   double FinalEner; ///<Minimum energy considered, in keV
00146   SIM_PARAM simpar; ///<Collection of parameters regarding to simulation
00147   }EXP_PARAM;
00148 
00149 /**Structure for passing info read during from input files*/
00150 typedef struct EXTRAINFO  
00151   {
00152    char *SampleFileNm;  ///<pointer to name of sample definition file
00153    char *FilterFileNm;  ///<pointer to name of filter definition file
00154    char *XYldFileNm;    ///<pointer to name of detector efficiency calibration file
00155    char *AreasFileNm;   ///<pointer to name of areas file (also used for calculation flags)
00156    char *DBpath;        ///<pointer to Path name for databases (for stoppings, FCK coefs, etc)
00157    char *OutputFileNm;   ///<pointer to name of output file
00158    int AreasFormat;     ///<Defines the format of the Areas File. 1 for raw areas format, 2 for DATTPIXE Areas Format
00159    int WantOutputfile;  ///<Flag indicating if an Output file is required (1) or not (0)
00160    int DoCalibration;   ///<Flag indicating that the input data is for performing a calibration (1) or not (0)
00161   }EXTRAINFO;
00162 
00163 /**Peak areas*/
00164 typedef struct PeakArea {
00165   ChemSymb symb;   ///< Symbol
00166   atomicnumber atnum; ///<Atomic number
00167   int area[5][3]; ///< Peak areas: [0][j]-->K  [1][j]-->LIII , [2][j]-->LII , [3][j]-->LI , [4][j]-->M
00168   int erro[5][3]; ///< Error in peak areas (same convention as in ".area")
00169 } PeakArea;
00170 
00171 /**L,K & M lines intensities*/
00172 typedef struct CalibYld {
00173   double K_[4];   ///< K Lines 
00174   double L_[4][4];   ///< 1-LIII , 2-LII , 3-LI 
00175   double M_[4];   ///< M  Lines
00176 } CalibYld;
00177 
00178 /** Emission energies and other ... 
00179 @todo document this structure*/
00180 typedef struct XrayYield {
00181   ChemSymb symb;
00182   atomicnumber atnum;
00183   CalibYld ener;
00184   CalibYld XYld;
00185 } XrayYield;
00186 
00187 /**Fluorescence and Coster-Kronig coefficients*/
00188 typedef struct FluorCKCoef {
00189   atomicnumber atnum;  
00190   double w[10];         
00191   double ck[3];        
00192   CalibYld k;
00193 } FluorCKCoef;
00194 
00195 /**@todo document this structure*/
00196 typedef struct SRType {
00197   CalibYld AbsFact;
00198   CalibYld Xeq;
00199   CalibYld SFCr;
00200   CalibYld DetLim;
00201   CalibYld ConcErr;
00202   CalibYld MDV;
00203 } SRType;
00204 
00205 
00206 /**@todo document this structure*/
00207 typedef struct SecResType {
00208   PeakArea Pk;
00209   SRType SR;
00210 } SecResType;
00211 
00212 
00213 
00214 /**@todo document this structure*/
00215 typedef struct ESxType {
00216   double ep;
00217   double stpp;
00218   double x;
00219 } ESxType;
00220  
00221 /**Basic sample foil description*/
00222 typedef struct foil {
00223   double thick; ///<Thickness, in 1e15at/cm2
00224   int nfoilelm; ///<Number of elements (redundant with foil.comp.nelem)
00225   COMPOUND comp; ///<compound definition
00226 } foil;
00227 
00228  
00229 /**Variables that depend on the layer*/  
00230 typedef struct 
00231   {
00232   int NumOfTrc;                ///<Number of "traces" (i.e. elements)
00233   int *TrcAtnum;               ///<Atomic numbers for each trace.
00234   int dimTrans;                ///<Dimension of Transmission vector
00235   CalibYld *Trans;             ///<Vector of calculated transmission for the lines of each element in this layer.
00236 //  PeakArea *AreasArray;        ///<Experimental Peak Areas. Not used for simulation. @todo check if it makes sense to maintain it here
00237 //  XrayYield *ResArray;        ///<@todo description?
00238   XrayYield *ResYldArray;     ///<Results (areas and energies of the calculated lines)
00239   SecResType *SecResArray;    ///<Same as ResYldArray but for secondary fluorescence results only
00240   int *TrcUse;          ///<Pointer to an array of flags controlling if calcs should be done for each given "trace" of this layer
00241   int *NeedSFC;         ///<Pointer to an array of flags controlling if sec. fluor. calcs should be done for each given element present in the SAMPLE
00242   ESxType *ESxArray;    ///<Array containing Energy, Stopping and position of each sublayer. 
00243   int FESxlen;          ///<Number of sublayers created for this layer (dimension of ESxArray)
00244   CalibYld *SSTrsArray; ///<Array containing Transmission coefs for each element in each sublayer
00245   CalibYld *SSYldArray; ///<Array containing X-ray production cross section for each element in each sublayer
00246   double ThickIn;       ///<In-going path length for the layer (and the curren geometry)
00247   double FoilInEner;    ///<Energy (keV) when entering the layer
00248   double FoilOutEner;    ///<Energy (keV) when exiting the layer
00249   double absolutepos;   ///<absolute in-going path length, in 1e15at/cm2, from sample surface. i.e. sumation of ThickIn for all layers below this one.
00250   const foil *pFoil;          ///<pointer to the foil associated to this layer (see foil structure)
00251   } LYR;
00252   
00253 
00254 /**X-ray Absortion coefficients*/
00255 typedef struct AbsCoef {
00256   atomicnumber atnum;     ///<Atomic number
00257   double coefstd[23];     ///<23 absorption coefficients, each one for a given standard energy tabulated in stdtable[] in TotAbsor() 
00258   double enr[10];          ///<Array transition energies (energies at which a sharp edge in absorption occurs): from 1 to 9: K LI LII LIII MI MII MIII MIV MV  (absorption E is decreasing)
00259   double coefenr[9][2];   ///<abs coefs for E just below (coefenr[i][0]) and just above (coefenr[i][1]) the transition energy enr[i].
00260 } AbsCoef;
00261 
00262 /**Definition of absorber filter. The filter can be defined in 2 alternative ways: 1.- Just like a you define a sample,  layer by layer (the transmission being calculated for each relevant line). 2.- Defining its transmission function by giving a series of coefficients that fit T(E) (this was the default way in DATPIXE).*/   
00263 typedef struct FILTER {
00264   int nlyr;         ///<number of layers defining the sample
00265   int MaxZ;         ///<maximum atomic number of elements present in sample
00266   int *FilterElems;  ///<Vector of flags indicating if a certain element is present in the filter (1) or not (0)
00267   double geomcorr;  ///<Geometrical correction to foil thickness (e.g. if filter is tilted: geomcorr=1/cos(tilt))
00268   foil *foil;       ///<definition of each layer
00269   int dimTrans;     ///<Dimension of Transmission vector
00270   CalibYld *Trans;  ///<Vector of calculated transmission for the lines of each element in the whole filter
00271   int changes;      ///<Flag indicating that Filter has changed (1) or not (0)
00272   int FilterType;   ///<Flag indicating type of filter used: 0->No filter 1->Foil-defined  2->Calibrated 
00273   absorber Calibflt;///<Structure to hold coeffs if alternative #2 is chosen. @todo (not supported yet)
00274 }FILTER;
00275 
00276 /**Calculated areas for a given element.*/
00277 typedef struct CPIXERESULTS{
00278   int atnum;          ///<atomic number of the element
00279   CalibYld simareas;  ///<K ,L & M yields 
00280   CalibYld err;       ///<K ,L & M yield uncertainties 
00281 }CPIXERESULTS;
00282 
00283 
00284 
00285 /**Each element which fluoresces knows which elements (and which lines) are the origin*/
00286 typedef struct SFCListElem {
00287   atomicnumber Za;  ///<Z of emissor of primary X rays
00288   atomicnumber Zb;  ///<Z of the element that fluoresces (secondary emissor)
00289   int epri;         ///<primary number identifying the emission group: 0=K, 1=LIII, 2=LII, 3=LI, 4=M
00290   int esec;         ///<secondary number identifying the emission line inside the group 
00291   int abs;          ///<Index of absorption energies. From 1 to 9: K LI LII LIII MI MII MIII MIV MV
00292   double Ea;        ///<Energy of the emission line (in keV)
00293   double sigmaphoto;///<photoemission cross section (@todo units?)
00294 } SFCListElem;
00295 
00296 /*Constants and Initializators*/  
00297 /**Limits for atomic numbers to be used
00298  @todo These values should be read instead of hardcoded. Maybe in SIMPARAM?
00299  */
00300 //static atomicnumber minK = 11, maxK = 54, minL = 30, maxL = 92, minM = 60;
00301 
00302 static const atomicnumber minK = 11, maxK = 54, minL = 30, minM = 60;
00303 
00304 static const CalibYld CYldNul = { { 0.0, 0.0, 0.0, 0.0 },
00305                             { { 0.0, 0.0, 0.0, 0.0 },
00306                               { 0.0, 0.0, 0.0, 0.0 },
00307                               { 0.0, 0.0, 0.0, 0.0 },
00308                               { 0.0, 0.0, 0.0, 0.0 }},
00309                             { 0.0, 0.0, 0.0, 0.0 }};
00310   
00311 /**Names for X-ray emmission lines*/
00312 
00313 //static LineNmType LineNm = {
00314 static char LineNm[6][4][11] = {
00315   { "--", "  --      ", "  --      ", "  --      " },
00316   { "--", "K_alpha1,2", "K_beta_1  ", "K_beta_2  " },
00317   { "--", "L_alpha1,2", "L_beta_2  ", "L_l       " },
00318   { "--", "L_beta_1  ", "L_gamma_1 ", "L_eta     " },
00319   { "--", "L_beta_3  ", "L_beta_4  ", "L_gamma_3 " },
00320   { "--", "M_alpha1,2", "M_beta    ", "M_gamma   " }
00321 };
00322 
00323   
00324   
00325 /*Function declarations*/
00326   
00327 void testlib(int **testarray,int *n);
00328 
00329 int symbol2Z(char *symbol);
00330 
00331 int Z2mass(int Z, double *mass, char option);
00332 
00333 int readCOMPOUND(FILE *f, int nelem, COMPOUND *c);
00334 
00335 int readINPUT(const char *InputFileNm, EXP_PARAM *pexppar, EXTRAINFO *pExtraInfo);
00336 
00337 int readEXP_PARAM(char *ExpParFileNm, EXP_PARAM *pexppar);
00338 
00339 int readSIM_PARAM(char *SampleParamFileNm, SIM_PARAM *simpar,char *IniMatFileNm, char *FilterFileNm, char *XYldFileNm);
00340 
00341 int readsample(char *SampleDefFileNm, int *MaxZ, int *NFoil, foil **Sample);
00342 
00343 int readFilter(const char *FilterDefFileNm, FILTER *Filter);
00344 
00345 int createPresentElems(int MaxZ, int NFoil,const foil *MatArray, int **PresentElems);
00346 
00347 int readCalcFlags(const char *CalcFlagsFileNm, const int *PresentElems, int MaxZinsample, int AreasFormat, CPIXERESULTS **CalcFlags);
00348 
00349 int readXYld(const char *XYldFileNm, const int *PresentElems, const CPIXERESULTS *CalcFlags, int MaxZinsample, double *CalEner, XrayYield **XYldArray);
00350 
00351 int readFCK(char *FCKCoefFileNm, int MaxZinsample, FluorCKCoef **FCKCoefArray);
00352 
00353 int readAbsCoef(char *TotAbsCoefFileNm, int MaxZinsample, const int *PresentElems, const FILTER *Filter, AbsCoef **TotAbsCoefArray);
00354 
00355 int createSPTs(const char *path, const EXP_PARAM *pexp, int MaxZinsample, const int *PresentElems, double step, SPT **SPTArray);
00356 
00357 double getSP(const COMPOUND *pcmp, const SPT *SPTArray, double E);
00358 
00359 int initlyrarray(const EXP_PARAM *pexp, const foil *MatArray, const XrayYield *XYldArray, const AbsCoef *TotAbsCoefArray, const SPT *SPTArray, const int *PresentElems, int NFoil, int *NFoilUsed, LYR **plyrarray);
00360 
00361 int initlyr(const EXP_PARAM *pexp, const XrayYield *XYldArray, const AbsCoef *TotAbsCoefArray, const int *PresentElems, LYR *plyr, double *MACoef);
00362 
00363 int FilterTrans(const EXP_PARAM *pexp, const XrayYield *XYldArray, const AbsCoef *TotAbsCoefArray, const int *PresentElems, FILTER *Filter);
00364 
00365 int Transmission(const XrayYield *XYldArray, const AbsCoef *TotAbsCoefArray, int Z, const foil *pFoil, double geomcorr, const CalibYld *pTransOld, CalibYld *pTrans);
00366 
00367 double TotAbsor(const AbsCoef *TotAbsCoefArray, const COMPOUND *cmp, double Xray);
00368 
00369 double TotAbsor_elemental(int iener, const AbsCoef *Absco, int Z, double Xray);
00370 
00371 int createsublyrs(const EXP_PARAM *pexp, const foil *pMat, const SPT *SPTArray, double MajAbsCoef, LYR *plyr);
00372 
00373 int SSThick(const EXP_PARAM *pexp, const foil *pMat, const SPT *SPTArray, double MajAbsCoef, double LayerThickness, ESxType ESxin,
00374         ESxType *ESxfin, int *pFpos, double *thick, ESxType **ESxA);
00375 
00376 int integrate_Simpson(const EXP_PARAM *pexp, const AbsCoef *TotAbsCoefArray, const FluorCKCoef *FCKCoefArray, const XrayYield *XYldArray, int NFoilUsed, const FILTER *Filter, LYR *plyrarray, CalibYld *XYldSums);
00377 
00378 int Xprod(const AbsCoef *AbsC, const FluorCKCoef *pFCK, atomicnumber Z1, atomicnumber Z2, double M2 , double ener, CalibYld *XYld);
00379 
00380 double PaulX(double ener, atomicnumber z);
00381 
00382 double PaulX_y(double MeV, atomicnumber z);
00383 
00384 double PaulX_lp(double x, long p);
00385 
00386 int ReisX(const AbsCoef *AbsC, const FluorCKCoef *pFCK, double ener, atomicnumber z, double M2, double *sigmaXL);
00387 
00388 double ReisX_gs(FwTipo T, int SS, double ksis);
00389 
00390 double ReisX_hs(FwTipo T, int SS, double ksih, double thet);
00391 
00392 double ReisX_En(double z, int niu);
00393 
00394 double ReisX_polisec(int ssind, double kz, double tz);
00395 
00396 double ReisX_g(int ss, atomicnumber Zg, double xi);
00397 
00398 void PenInteg(atomicnumber atnumb, const CalibYld *absc, const ESxType *ESA,
00399         const CalibYld *YldA, const CalibYld *TrsA, const CalibYld *pTrs0, 
00400         int FExlen, int NeedSFC, int AllowXEqCalc,
00401         double x0, double CosInc,
00402         CalibYld *XYld, CalibYld *XSFCr, CalibYld *XYldxmed);
00403 
00404 double Simps(double a, double b, double fa, double fi, double fb);
00405 
00406 void deNormalize(const EXP_PARAM *pexp, const AbsCoef *AbsC, 
00407                const FluorCKCoef *pFCK, atomicnumber Z2, double M2,  double attfraction, 
00408                const CalibYld *pXYld, XrayYield *ResY, CalibYld *XYldSum);
00409 
00410 void freeFilter(FILTER *Filter);    
00411                
00412 void freeReusable(int NFoils, LYR **plyrarray, foil **MatArray, CalibYld **XYldSums );
00413 
00414 void safefree(void **ptr);
00415 
00416 void fprintCALIBYLD(FILE *f,const CalibYld *pCYld);
00417 
00418 int createSFCList(const EXP_PARAM *pexp, const int *PresentElems, const XrayYield *XYldArray, const AbsCoef *TotAbsCoefArray, int *dimList, SFCListElem **SFCList);
00419 
00420 int needSFC(atomicnumber Za, atomicnumber Zb, const CalibYld *enerA, const AbsCoef *AbsC, int *dimList, SFCListElem **SFCList );
00421 
00422 double sigmaphoto(int iabs, const AbsCoef *AbsC, atomicnumber Zb, double Xray);
00423 
00424 int inrange(double value,double min, double max);
00425 

Generated on Fri Jul 15 20:43:49 2005 for LibCPIXE API by  doxygen 1.3.9.1