/home66/gary/public_html/cloudy/c08_branch/source/atmdat_chianti.cpp

Go to the documentation of this file.

/* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and
* others.  For conditions of distribution and use see copyright notice in license.txt */
#include "cddefines.h"
#include "lines_service.h"
#include "physconst.h"
#include "taulines.h"
#include "trace.h"
#include "string.h"
#include "thirdparty.h"

#define DEBUGSTATE false

void atmdat_Chianti_readin( long intNS )
{
        DEBUG_ENTRY( "atmdat_Chianti_readin()" );

        int intCS,intCollIndex = -10000,intsplinepts,intTranType,intxs;
        /* type is set to 0 for non chianti and 1 for chianti*/
        long int i,j,nMolLevs,intCollTran,intcollindex,ipLo,ipHi;
        FILE *atmolLevDATA , *atmolTraDATA=NULL,*atmolEleColDATA=NULL,*atmolProColDATA=NULL;
        realnum  fstatwt,fenergyK,fenergyWN,fWLAng,fenergy,feinsteina;
        double fScalingParam,fEnergyDiff,fGF,*xs,*spl,*spl2;

        //Energy Column Start
        //Statistical Weight Start
        //Lower Level Start
        //Upper Level Start
        //Einstein A Start
        //Wavelength of transition start
        long ECS,SWS,LLS,ULS,EAS,WOTS;

        char chLine[FILENAME_PATH_LENGTH_2] ,
                chEnFilename[FILENAME_PATH_LENGTH_2],
                chTraFilename[FILENAME_PATH_LENGTH_2],
                chEleColFilename[FILENAME_PATH_LENGTH_2],
                chProColFilename[FILENAME_PATH_LENGTH_2];

        char *chDesired;
        realnum *atmolStatesEnergy;
        long int *intNewIndex=NULL,*intOldIndex=NULL;
        int interror;
        bool lgProtonData=false,lgEneLevOrder;

        /*This is to identify where the fields start */
        /*In Chianti(1),the energy field in cm-1 starts at 40 and has a width 15*/
        ECS = 71;/*We use the theoretical energy in cm-1*/
        SWS = 38;
        LLS = 1;
        ULS = 6;
        EAS = 41;
        WOTS = 11;

        /*For the CHIANTI DATABASE*/
        /*First Opening the energy levels file*/
        strcpy( chEnFilename , Species[intNS].chptrSpName );                    
        strcat( chEnFilename , ".elvlc");
        uncaps( chEnFilename );
        if(DEBUGSTATE)
                printf("The data file is %s \n",chEnFilename);

        /*Open the files*/
        if( trace.lgTrace )
                fprintf( ioQQQ," moldat_readin opening %s:",chEnFilename);

        atmolLevDATA = open_data( chEnFilename, "r" );

        /*Second:Opening the transition probabilities file*/
        strcpy( chTraFilename , Species[intNS].chptrSpName );                   
        strcat( chTraFilename , ".wgfa");
        uncaps( chTraFilename );
        if(DEBUGSTATE)
                printf("The data file is %s \n",chTraFilename);

        /*Open the files*/
        if( trace.lgTrace )
                fprintf( ioQQQ," moldat_readin opening %s:",chTraFilename);

        atmolTraDATA = open_data( chTraFilename, "r" );

        /*Third: Opening the electron collision data*/
        strcpy( chEleColFilename , Species[intNS].chptrSpName );                        
        strcat( chEleColFilename , ".splups");
        uncaps( chEleColFilename );
        if(DEBUGSTATE)
                printf("The data file is %s \n",chEleColFilename);
        /*Open the files*/
        if( trace.lgTrace )
                fprintf( ioQQQ," moldat_readin opening %s:",chEleColFilename);

        atmolEleColDATA = open_data( chEleColFilename, "r" );

        /*Fourth: Opening the proton collision data*/
        strcpy( chProColFilename , Species[intNS].chptrSpName );                        
        strcat( chProColFilename , ".psplups");
        uncaps( chProColFilename );
        if(DEBUGSTATE)
                printf("The data file is %s \n",chProColFilename);
        /*Open the files*/
        if( trace.lgTrace )
                fprintf( ioQQQ," moldat_readin opening %s:",chProColFilename);
        /*We will set a flag here to indicate if the proton collision strengths are available */
        if( ( atmolProColDATA = fopen( chProColFilename , "r" ) ) != NULL )
        {
                lgProtonData = true;
                fclose( atmolProColDATA );
                atmolProColDATA = NULL;
        }
        else
        {
                lgProtonData = false;
        }

        nMolLevs = 0;
        lgEneLevOrder = true;
        while( read_whole_line( chLine , (int)sizeof(chLine) , atmolLevDATA ) != NULL )
        {       
                chDesired = strtok(chLine," \t");
                intCS = atoi(chDesired);
                if (intCS == -1)
                {       
                        break;
                }
                else
                        nMolLevs++;
        }

        if(DEBUGSTATE)
                printf("The number of energy levels is %li",nMolLevs);

        if( nMolLevs <= 0 )
        {
                fprintf( ioQQQ, "The number of energy levels is non-positive in datafile %s.\n", chEnFilename );
                fprintf( ioQQQ, "The file must be corrupted.\n" );
                cdEXIT( EXIT_FAILURE );
        }

        Species[intNS].numLevels_max = nMolLevs;
        Species[intNS].numLevels_local = nMolLevs;
        /*Malloc the energy array to check if the energy levels are in order*/
        atmolStatesEnergy = (realnum*)MALLOC((unsigned long)(nMolLevs)*sizeof(realnum));
        /*malloc the States array*/
        atmolStates[intNS] = (quantumState *)MALLOC((size_t)(nMolLevs)*sizeof(quantumState));
        /*malloc the Transition array*/
        atmolTrans[intNS] = (transition **)MALLOC(sizeof(transition*)*(unsigned long)nMolLevs);
        atmolTrans[intNS][0] = NULL;
        for( ipHi = 1; ipHi < nMolLevs; ipHi++)
        {
                atmolTrans[intNS][ipHi] = (transition *)MALLOC(sizeof(transition)*(unsigned long)ipHi);
                for(ipLo = 0; ipLo < ipHi; ipLo++)
                {
                        atmolTrans[intNS][ipHi][ipLo].Lo = &atmolStates[intNS][ipLo];
                        atmolTrans[intNS][ipHi][ipLo].Hi = &atmolStates[intNS][ipHi];
                        atmolTrans[intNS][ipHi][ipLo].Emis = &DummyEmis;
                }
        }

        for( intcollindex = 0; intcollindex<NUM_COLLIDERS; intcollindex++ )
        {
                CollRatesArray[intNS][intcollindex] = NULL;
        }
        /*Allocating space just for the electron*/
        CollRatesArray[intNS][0] = (double**)MALLOC((nMolLevs)*sizeof(double*));
        for( ipHi = 0; ipHi<nMolLevs; ipHi++ )
        {
                CollRatesArray[intNS][0][ipHi] = (double*)MALLOC((unsigned long)(nMolLevs)*sizeof(double));
                for( ipLo = 0; ipLo<nMolLevs; ipLo++)
                {
                        CollRatesArray[intNS][0][ipHi][ipLo] = 0.0;
                }
        }
        /*Allocating space for the proton*/
        if(lgProtonData)
        {
                CollRatesArray[intNS][1] = (double**)MALLOC((nMolLevs)*sizeof(double*));
                for( ipHi = 0; ipHi<nMolLevs; ipHi++ )
                {
                        CollRatesArray[intNS][1][ipHi] = (double*)MALLOC((unsigned long)(nMolLevs)*sizeof(double));
                        for( ipLo = 0; ipLo<nMolLevs; ipLo++)
                        {
                                CollRatesArray[intNS][1][ipHi][ipLo] = 0.0;
                        }
                }
        }
        /*Rewind the energy levels files*/
        if( fseek( atmolLevDATA , 0 , SEEK_SET ) != 0 )
        {
                fprintf( ioQQQ, " moldat_readin could not rewind %s.\n",chEnFilename);
                cdEXIT(EXIT_FAILURE);
        }


        /*Reading in the energy and checking that they are in order*/
        for( i=0; i<nMolLevs; i++ )
        {
                if(read_whole_line( chLine , (int)sizeof(chLine)  , atmolLevDATA ) == NULL )
                {
                        fprintf( ioQQQ, " The data file %s is corrupted .\n",chEnFilename);
                        cdEXIT(EXIT_FAILURE);
                }
                fenergy = (realnum) atof(&chLine[ECS]);
                atmolStatesEnergy[i] = fenergy;

                /*To check if energy levels are in order*/
                /*If the energy levels are not in order a flag is set*/
                if(i>0)
                {
                        if (atmolStatesEnergy[i] < atmolStatesEnergy[i-1])
                        {
                                lgEneLevOrder = false;  
                        }
                }
        }

        if(DEBUGSTATE)
        {
                for(i=0;i<nMolLevs;i++)
                {
                        printf("The %ld value of the unsorted energy array is %f \n",i,atmolStatesEnergy[i]);
                }
        }

        intNewIndex = (long int*)MALLOC((unsigned long)(nMolLevs)* sizeof(long int));

        if(lgEneLevOrder == false)
        {
                /*If the energy levels are not in order and the flag is set(lgEneLevOrder=FALSE)
                then we sort the energy levels to find the relation matrix between the old and new indices*/
                intOldIndex = (long int*)MALLOC((unsigned long)(nMolLevs)* sizeof(long int));
                /*We now do a modified quick sort*/
                spsort(atmolStatesEnergy,nMolLevs,intOldIndex,2,&interror);
                /*intNewIndex has the key to map*/
                for( i=0; i<nMolLevs; i++ )
                {
                        ASSERT( intOldIndex[i] < nMolLevs );
                        intNewIndex[intOldIndex[i]] = i;        
                }
                if(DEBUGSTATE)
                {
                        for(i=0;i<nMolLevs;i++)
                        {
                                printf("The %ld value of the relation matrix is %ld \n",i,intNewIndex[i]);
                        }
                        for( i=0; i<nMolLevs; i++ )
                        {
                                printf("The %ld value of the sorted energy array is %f \n",i,atmolStatesEnergy[i]);
                        }
                }

                free( intOldIndex );
        }
        else
        {
                /* if energies already in correct order, new index is same as original. */
                for( i=0; i<nMolLevs; i++ )
                {
                        intNewIndex[i] = i;     
                }
        }

        /* insist that there the intNewIndex values are unique */
        for( i=0; i<nMolLevs; i++ )
        {
                for( j=i+1; j<nMolLevs; j++ )
                {
                        ASSERT( intNewIndex[i] != intNewIndex[j] );     
                }
        }


        /*After we read in the energies we rewind the energy levels file */
        if( fseek( atmolLevDATA , 0 , SEEK_SET ) != 0 )
        {
                fprintf( ioQQQ, " moldat_readin could not rewind %s.\n",chEnFilename);
                cdEXIT(EXIT_FAILURE);
        }

        for( i=0; i<nMolLevs; i++ )
        {
                long j = intNewIndex[i];

                if(read_whole_line( chLine , (int)sizeof(chLine)  , atmolLevDATA ) == NULL )
                {
                        fprintf( ioQQQ, " The data file %s is corrupted .\n",chEnFilename);
                        cdEXIT(EXIT_FAILURE);
                }
                /*information needed for label*/
                strcpy(atmolStates[intNS][j].chLabel, "    ");
                strncpy(atmolStates[intNS][j].chLabel,Species[intNS].chptrSpName, 4);
                atmolStates[intNS][j].chLabel[4] = '\0';

                fstatwt = (realnum)atof(&chLine[SWS]);
                if(fstatwt <= 0.)
                {
                        fprintf( ioQQQ, " A positive non zero value is expected for the statistical weight .\n");
                        fprintf( ioQQQ, " The data file %s is corrupted .\n",chEnFilename);
                        cdEXIT(EXIT_FAILURE);

                }
                atmolStates[intNS][j].g = fstatwt;
                fenergy = (realnum) atof(&chLine[ECS]);
                atmolStates[intNS][j].energy = fenergy;

                if(DEBUGSTATE)
                {
                        fprintf( ioQQQ, "The %ld value of g after populating is %f \n",
                                j,atmolStates[intNS][j].g);
                        fprintf( ioQQQ, "The %ld value of energy after populating is %f \n",
                                j,atmolStates[intNS][j].energy);
                }
        }


        /* fill in all transition energies, can later overwrite for specific radiative transitions */
        for( i=1; i<nMolLevs; i++ )
        {
                for( j=0; j<i; j++ )
                {
                        fenergyWN = MAX2( (realnum)1E-10, atmolStates[intNS][i].energy - atmolStates[intNS][j].energy );
                        fenergyK = (realnum)(fenergyWN*T1CM);

                        atmolTrans[intNS][i][j].EnergyK = fenergyK;
                        atmolTrans[intNS][i][j].EnergyWN = fenergyWN;
                        atmolTrans[intNS][i][j].EnergyErg = (realnum)ERG1CM *fenergyWN;
                        atmolTrans[intNS][i][j].WLAng = (realnum)(1e+8/fenergyWN/RefIndex(fenergyWN));
                }
        }

        /************************************************************************/
        /*** Read in the level numbers, Einstein A and transition wavelength  ***/
        /************************************************************************/
        if(read_whole_line( chLine , (int)sizeof(chLine)  ,atmolTraDATA ) == NULL )
        {
                fprintf( ioQQQ, " The data file %s is corrupted .\n",chTraFilename);
                cdEXIT(EXIT_FAILURE);
        }

        while(atoi(chLine)!= -1)
        {
                ipLo = atoi(&chLine[LLS])-1;
                ipHi = atoi(&chLine[ULS])-1;

                /* translate to properly sorted indices */
                ipLo = intNewIndex[ipLo];
                ipHi = intNewIndex[ipHi];

                ASSERT( ipHi > ipLo );

                fWLAng = (realnum)atof(&chLine[WOTS]);

                /* \todo 2 Chianti labels the H 1 2-photon transition as z wavelength of zero.  
                 * Should we just ignore all of the wavelengths in this file and use the
                 * difference of level energies instead. */
                if( fWLAng == 0. )
                {
                        fWLAng = (realnum)1e8/
                                (atmolStates[intNS][ipHi].energy - atmolStates[intNS][ipLo].energy);
                }

                feinsteina = (realnum)atof(&chLine[EAS]);
                atmolTrans[intNS][ipHi][ipLo].Emis 
                                = AddLine2Stack(feinsteina, &atmolTrans[intNS][ipHi][ipLo]);
                atmolTrans[intNS][ipHi][ipLo].WLAng = fWLAng;
                fenergyWN = (realnum)(1e+8/fWLAng);
                fenergyK = (realnum)(fenergyWN*T1CM);

                atmolTrans[intNS][ipHi][ipLo].EnergyK = fenergyK;
                atmolTrans[intNS][ipHi][ipLo].EnergyWN = fenergyWN;
                atmolTrans[intNS][ipHi][ipLo].EnergyErg = (realnum)ERG1CM *fenergyWN;
                atmolTrans[intNS][ipHi][ipLo].WLAng = (realnum)(1e+8/fenergyWN/RefIndex(fenergyWN));

                ASSERT( !isnan( atmolTrans[intNS][ipHi][ipLo].EnergyK ) );

                if(DEBUGSTATE)
                {
                        fprintf( ioQQQ, "The lower level is %ld \n",ipLo);
                        fprintf( ioQQQ, "The upper level is %ld \n",ipHi);
                        fprintf( ioQQQ, "The Einstein A is %f \n",atmolTrans[intNS][ipHi][ipLo].Emis->Aul);
                        fprintf( ioQQQ, "The wavelength of the transition is %f \n",atmolTrans[intNS][ipHi][ipLo].WLAng );
                }

                if(read_whole_line( chLine , (int)sizeof(chLine) , atmolTraDATA ) == NULL )
                {
                        fprintf( ioQQQ, " Data is expected on this line.\n");
                        fprintf( ioQQQ, " The data file %s is corrupted .\n",chTraFilename);
                        cdEXIT(EXIT_FAILURE);
                }
        }

        /* Malloc space for splines */

        AtmolCollSplines[intNS] = (CollSplinesArray***)MALLOC(nMolLevs *sizeof(CollSplinesArray**));
        for( i=0; i<nMolLevs; i++ )
        {
                AtmolCollSplines[intNS][i] = (CollSplinesArray **)MALLOC((unsigned long)(nMolLevs)*sizeof(CollSplinesArray *));
                for( j=0; j<nMolLevs; j++ )
                {
                        AtmolCollSplines[intNS][i][j] = 
                                (CollSplinesArray *)MALLOC((unsigned long)(NUM_COLLIDERS)*sizeof(CollSplinesArray ));

                        for( long k=0; k<NUM_COLLIDERS; k++ )
                        {
                                /* initialize all spline variables */
                                AtmolCollSplines[intNS][i][j][k].collspline = NULL;
                                AtmolCollSplines[intNS][i][j][k].SplineSecDer = NULL;
                                AtmolCollSplines[intNS][i][j][k].nSplinePts = -1; 
                                AtmolCollSplines[intNS][i][j][k].intTranType = -1;
                                AtmolCollSplines[intNS][i][j][k].gf = BIGDOUBLE;
                                AtmolCollSplines[intNS][i][j][k].EnergyDiff = BIGDOUBLE;
                                AtmolCollSplines[intNS][i][j][k].ScalingParam = BIGDOUBLE;
                        }
                }
        }

        fixit(); /* This is where the loop of colliders should begin */

        /*********************************************/
        /*** Read in the electron collisional data ***/
        /*********************************************/
        intCollIndex = 0;
        intCollTran = 0;
        if(read_whole_line( chLine , (int)sizeof(chLine)  ,atmolEleColDATA ) == NULL )
        {
                fprintf( ioQQQ, " The data file %s is corrupted .\n",chEleColFilename);
                cdEXIT(EXIT_FAILURE);
        }
        /*We need to find the number of collisional transitions */
        while(atoi(chLine)!= -1)
        {
                intCollTran++;
                if(read_whole_line( chLine , (int)sizeof(chLine) , atmolEleColDATA ) == NULL )
                {
                        fprintf( ioQQQ, " Data is expected on this line.\n");
                        fprintf( ioQQQ, " The data file %s is corrupted .\n",chEleColFilename);
                        cdEXIT(EXIT_FAILURE);

                }
        }

        if(DEBUGSTATE)
                printf("\n The number of collisional transitions is %li",intCollTran);

        /*We rewind the file*/
        if( fseek( atmolEleColDATA , 0 , SEEK_SET ) != 0 )
        {
                fprintf( ioQQQ, " moldat_readin could not rewind %s.\n",chEleColFilename);
                cdEXIT(EXIT_FAILURE);
        }
        /*Dummy string used for convenience*/
        strcpy(chLine,"A");

        if(read_whole_line( chLine , (int)sizeof(chLine)  ,atmolEleColDATA ) == NULL )
        {
                fprintf( ioQQQ, " The data file %s is corrupted .\n",chEleColFilename);
                cdEXIT(EXIT_FAILURE);
        }

        /*If the file exists,the first line cannot be "-1"*/
        if(atoi(chLine)== -1)
        {
                fprintf( ioQQQ, " If the file %s exists,the first line cannot be -1 .\n",chEleColFilename);
                cdEXIT(EXIT_FAILURE);
        }

        while(atoi(chLine)!= -1)
        {
                i = 1;
                bool lgEOL;

                (void)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL ); // intAtomicNo
                (void)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL ); // intIonStage
                /* level indices */
                ipLo = (long)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL ) - 1;
                ipHi = (long)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL ) - 1;
                intTranType = (long)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL );
                /*gf*/
                fGF = (realnum)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL );
                /*Energy Difference*/
                fEnergyDiff = (realnum)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL );
                /*Scaling Parameter*/
                fScalingParam = (realnum)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL );

                /* translate to properly sorted indices */
                ipLo = intNewIndex[ipLo];
                ipHi = intNewIndex[ipHi];

                ASSERT( ipLo < nMolLevs );
                ASSERT( ipHi < nMolLevs );
                ASSERT( AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].collspline == NULL );
                ASSERT( AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].SplineSecDer == NULL );

                fixit(); /* use a macro for the 5 and 9 that appear here. */
                /*We malloc the space here*/
                AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].collspline = 
                        (double *)MALLOC((unsigned long)(9)*sizeof(double));
                AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].SplineSecDer = 
                        (double *)MALLOC((unsigned long)(9)*sizeof(double));

                /* always read at least 5 */
                for( intsplinepts=0; intsplinepts<10; intsplinepts++ )
                {
                        double temp = (double)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL );

                        if( lgEOL )
                        {
                                /* there should be 5 or 9 spline points.  If we got EOL, 
                                 * insist that we were trying to read the 6th or 10th. */
                                if( (intsplinepts != 5) && (intsplinepts != 9) )
                                {
                                        fprintf( ioQQQ, " There should have been 5 or 9 spline points.  Sorry.\n" );
                                        fprintf( ioQQQ, "string==%s==\n" ,chLine );
                                        cdEXIT(EXIT_FAILURE);
                                }
                                else
                                        break;
                        }
                        else 
                        {
                                ASSERT( intsplinepts < 9 );
                                AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].collspline[intsplinepts] = temp;
                        }
                }

                ASSERT( (intsplinepts == 5) || (intsplinepts == 9) );

                /*The zeroth element contains the number of spline points*/             
                AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].nSplinePts = intsplinepts; 
                /*Transition type*/
                AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].intTranType = intTranType;
                /*gf value*/
                AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].gf = fGF;
                /*Energy difference between two levels in Rydbergs*/
                AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].EnergyDiff = fEnergyDiff;
                /*Scaling parameter C*/
                AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].ScalingParam = fScalingParam;

                /*Once the spline points have been filled,fill the second derivatives structure*/
                /*Creating spline points array*/
                xs = (double *)MALLOC((unsigned long)intsplinepts*sizeof(double));
                spl = (double *)MALLOC((unsigned long)intsplinepts*sizeof(double));
                spl2 = (double *)MALLOC((unsigned long)intsplinepts*sizeof(double));

                if(intsplinepts == 5)
                {
                        for(intxs=0;intxs<5;intxs++)
                        {
                                xs[intxs] = 0.25*intxs;
                                spl[intxs] = AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].collspline[intxs];
                        }
                }
                else if(intsplinepts == 9)
                {
                        for(intxs=0;intxs<9;intxs++)
                        {
                                xs[intxs] = 0.125*intxs;
                                spl[intxs] = AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].collspline[intxs];       
                        }
                }
                else
                        TotalInsanity();

                spline(xs, spl,intsplinepts,2e31,2e31,spl2);

                /*Filling the second derivative structure*/
                for(i=0;i<intsplinepts;i++)
                {
                        AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].SplineSecDer[i] = spl2[i];
                }

                free( xs );
                free( spl );
                free( spl2 );

                if(read_whole_line( chLine , (int)sizeof(chLine)  ,atmolEleColDATA ) == NULL )
                {
                        fprintf( ioQQQ, " The data file %s is corrupted .\n",chEleColFilename);
                        cdEXIT(EXIT_FAILURE);
                }
        }

        /*****************************************/
        /*** Print the AtmolCollSplines Values ***/
        /*****************************************/

        if(DEBUGSTATE)
        {
                intNS = 0;
                intCollIndex = 0;
                ipHi = 1;
                ipLo = 0;

                fprintf( ioQQQ, "Collisional Data: %li\t%li\t%f\t%f\t%f",
                        AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].nSplinePts,
                        AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].intTranType,
                        AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].gf,
                        AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].EnergyDiff,
                        AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].ScalingParam );
                for( j=0; j< AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].nSplinePts; j++)
                {
                        fprintf( ioQQQ, "\t%f",AtmolCollSplines[intNS][ipHi][ipLo][intCollIndex].collspline[j]);
                }
                fprintf( ioQQQ, "\n" );
        }

        free( atmolStatesEnergy );
        free( intNewIndex );

        fclose( atmolLevDATA );
        fclose( atmolTraDATA );
        fclose( atmolEleColDATA );

        return;
}

---