cloudy trunk
Loading...
Searching...
No Matches
save_average.cpp
Go to the documentation of this file.
1/* This file is part of Cloudy and is copyright (C)1978-2013 by Gary J. Ferland and
2 * others. For conditions of distribution and use see copyright notice in license.txt */
3/*save_average routine to read in list of species to output as averages */
4#include "cddefines.h"
5#include "cddrive.h"
6#include "input.h"
7#include "elementnames.h"
8#include "save.h"
9#include "parser.h"
10
11/* return value is number of lines, -1 if file could not be opened */
13 Parser &p,
14 /* the file we will write to */
15 long int ipPun,
16 char *chHeader)
17{
18 long int i;
19 long int nLine;
20
21 DEBUG_ENTRY( "parse_save_average()" );
22
23 char chCap[INPUT_LINE_LENGTH];
24 char chTemp[INPUT_LINE_LENGTH];
25
26 /* this is index for first line we will read. use this to count
27 * total number of species */
28 nLine = input.nRead;
29 /* very first time this routine is called, chJob is "READ" and we read
30 * in lines from the input stream. The species labels and other information
31 * are store in the save structure. These are output in a later call
32 * to this routine with argument PUNCH */
33
34 /* number of lines we will save */
35 save.nAverageList[ipPun] = 0;
36
37 /* get the next line, and check for eof */
38 p.getline();
39 if( p.m_lgEOF )
40 {
41 fprintf( ioQQQ,
42 " Punch average hit EOF while reading list; use END to end list.\n" );
44 }
45
46 /* keep reading until we hit END */
47 while( p.strcmp( "END" ) != 0 )
48 {
49 /* count number of species we will save */
50 ++save.nAverageList[ipPun];
51
52 /* get next line and check for eof */
53 p.getline();
54 if( p.m_lgEOF )
55 {
56 fprintf( ioQQQ, " save averages hit EOF while reading species list; use END to end list.\n" );
58 }
59 }
60/*# define PADEBUG*/
61# ifdef PADEBUG
62 fprintf(ioQQQ , "DEBUG save_average %li species read in.\n",
63 save.nAverageList[ipPun] );
64# endif
65
66 /* now create space that will be needed to hold these arrays */
67
68 if( save.ipPnunit[ipPun] == NULL )
69 {
70 /* make sure the memory from a previous call is freed */
71 save.SaveAverageFree(ipPun);
72
73 /* nAverageIonList is set of ions for averages */
74 save.nAverageIonList[ipPun].resize(save.nAverageList[ipPun]);
75
76 /* nAverage2ndPar is set of second parameters for averages */
77 save.nAverage2ndPar[ipPun].resize(save.nAverageList[ipPun]);
78
79 /* chAverageType is label for type of average */
80 save.chAverageType[ipPun].resize(save.nAverageList[ipPun]);
81
82 /* chAverageSpeciesLabel is label for species */
83 save.chAverageSpeciesLabel[ipPun].resize(save.nAverageList[ipPun]);
84
85 for( i=0; i < save.nAverageList[ipPun]; ++i )
86 {
87 /* create space for labels themselves */
88 save.chAverageType[ipPun][i] = new char[5];
89 /* chAverageSpeciesLabel is label for species */
90 save.chAverageSpeciesLabel[ipPun][i] = new char[5];
91 }
92 }
93
94 /* reset array input read to first of the species lines */
95 // Evil input rewind CodeSmell -- grrw
96 input.nRead = nLine;
97
98# ifdef PADEBUG
99 fprintf(ioQQQ , "DEBUG save_average %li species read in.\n",
100 save.nAverageList[ipPun] );
101# endif
102
103 /* reread the input lines and save the data */
104 p.getline();
105 if( p.m_lgEOF )
106 {
107 fprintf( ioQQQ,
108 " Punch average hit EOF while reading list; use END to end list.\n" );
110 }
111
112 /* use this to count number of species, and will assert equal to
113 * total malloced above */
114 nLine = 0;
115 /* keep reading until we hit END */
116 while( p.strcmp("END" ) != 0 )
117 {
118 /* count number of species we will save */
119 ++nLine;
120 if( p.nMatch("TEMP" ))
121 {
122 /* temperature */
123 strcpy( save.chAverageType[ipPun][nLine-1] , "TEMP" );
124 }
125 else if( p.nMatch("COLU" ))
126 {
127 /* column density */
128 strcpy( save.chAverageType[ipPun][nLine-1] , "COLU" );
129 }
130 else if( p.nMatch("IONI" ))
131 {
132 /* ionization fraction */
133 strcpy( save.chAverageType[ipPun][nLine-1] , "IONI" );
134 }
135 else
136 {
137 fprintf(ioQQQ,"PROBLEM one of the jobs TEMPerature, COLUmn density, or IONIzation, must appear.\n");
139 }
140
141 /* get element name, a string we shall pass on to the routine
142 * that computes final quantities */
143 if( (i = p.GetElem( )) < 0 )
144 {
145 /* no name found */
146 fprintf(ioQQQ, "save average did not an element on this line, sorry %s\n",
147 chCap );
149 }
150 strcpy( save.chAverageSpeciesLabel[ipPun][nLine-1] , elementnames.chElementNameShort[i]);
151
152 /* now get ionization stage */
153 save.nAverageIonList[ipPun][nLine-1] = (int) p.FFmtRead();
154 if( p.lgEOL() )
155 {
156 /* error - needed that ionization stage */
157 p.NoNumb("ionization stage" );
158 }
159
160 /* look for volume keyword, otherwise will be radius
161 * only used for some options */
162 if( p.nMatch( "VOLU" ) )
163 {
164 /* volume */
165 save.nAverage2ndPar[ipPun][nLine-1] = 1;
166 }
167 else
168 {
169 /* radius */
170 save.nAverage2ndPar[ipPun][nLine-1] = 0;
171 }
172
173 /* get next line and check for eof */
174 p.getline();
175 if( p.m_lgEOF )
176 {
177 fprintf( ioQQQ, " save averages hit EOF while reading species list; use END to end list.\n" );
179 }
180 }
181
182 /* these must equal or we have a major logical error */
183 ASSERT( nLine == save.nAverageList[ipPun]);
184
185# ifdef PADEBUG
186 for( i=0; i<nLine ; ++i )
187 {
188 fprintf(ioQQQ, "PDDEBUG %s %s %i %i\n",
189 save.chAverageType[ipPun][i],
190 save.chAverageSpeciesLabel[ipPun][i] ,
191 save.nAverageIonList[ipPun][i] ,
192 save.nAverage2ndPar[ipPun][i] );
193 }
194# endif
195
196 /* save headers */
197 sprintf(chHeader, "#averages");
198 for( i=0; i<nLine ; ++i )
199 {
200 sprintf(chTemp, "\t %s %s %i %i",
201 save.chAverageType[ipPun][i],
202 save.chAverageSpeciesLabel[ipPun][i] ,
203 save.nAverageIonList[ipPun][i] ,
204 save.nAverage2ndPar[ipPun][i] );
205 strcat( chHeader, chTemp );
206 }
207 strcat( chHeader, "\n");
208}
209
211 /* the file we will write to */
212 long int ipPun)
213{
214 long int i;
215
216 DEBUG_ENTRY( "save_average()" );
217
218 /* do the output */
219 for( i=0; i<save.nAverageList[ipPun] ; ++i )
220 {
221 double result;
222 char chWeight[7];
223 if( save.nAverage2ndPar[ipPun][i] == 0 )
224 strcpy( chWeight , "RADIUS");
225 else
226 strcpy( chWeight , "VOLUME");
227
228 if( strncmp( save.chAverageType[ipPun][i] , "TEMP" , 4 ) == 0)
229 {
230 /* temperature */
231 if( cdTemp(
232 save.chAverageSpeciesLabel[ipPun][i] ,
233 save.nAverageIonList[ipPun][i] ,
234 &result ,
235 chWeight ) )
236 {
237 fprintf( ioQQQ, " save average temperature could not identify the species.\n" );
239 }
240 /* will report log of temperature */
241 result = log10( result );
242 }
243 else if( strncmp( save.chAverageType[ipPun][i] , "IONI" , 4 ) == 0)
244 {
245 /* ionization fraction
246 * H2 is a special case, HYDRO 0 requests
247 * the H2 fraction, n(H2)/n(H) */
248 if( strncmp( "HYDR" ,
249 save.chAverageSpeciesLabel[ipPun][i] , 4)==0 &&
250 save.nAverageIonList[ipPun][i]== 0 )
251 strncpy( save.chAverageSpeciesLabel[ipPun][i],
252 "H2 " , 4 );
253 if( cdIonFrac(
254 save.chAverageSpeciesLabel[ipPun][i] ,
255 save.nAverageIonList[ipPun][i] ,
256 &result ,
257 chWeight ,
258 false
259 ) )
260 {
261 fprintf( ioQQQ, " save average ionization fraction could not identify the species.\n" );
263 }
264 /* will report log of ionization fraction */
265 result = log10( result );
266 }
267 else if( strncmp( save.chAverageType[ipPun][i] , "COLU" , 4 ) == 0)
268 {
269 /* column density */
270 if( cdColm(
271 save.chAverageSpeciesLabel[ipPun][i] ,
272 save.nAverageIonList[ipPun][i] ,
273 &result ) )
274 {
275 fprintf( ioQQQ, " save average column density fraction could not identify the species.\n" );
277 }
278 /* will report log of column density */
279 result = log10( result );
280 }
281 else
283
284 fprintf(save.ipPnunit[ipPun], "\t %.3f", result );
285 }
286 fprintf(save.ipPnunit[ipPun], "\n");
287}
FILE * ioQQQ
Definition cddefines.cpp:7
#define ASSERT(exp)
Definition cddefines.h:578
const int INPUT_LINE_LENGTH
Definition cddefines.h:254
#define EXIT_FAILURE
Definition cddefines.h:140
#define cdEXIT(FAIL)
Definition cddefines.h:434
NORETURN void TotalInsanity(void)
Definition service.cpp:886
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
int cdColm(const char *chLabel, long int ion, double *theocl)
Definition cddrive.cpp:636
int cdIonFrac(const char *chLabel, long int IonStage, double *fracin, const char *chWeight, bool lgDensity)
Definition cddrive.cpp:1072
int cdTemp(const char *chLabel, long int IonStage, double *TeMean, const char *chWeight)
Definition cddrive.cpp:1602
bool getline(void)
Definition parser.cpp:164
long int GetElem(void) const
Definition parser.cpp:209
double FFmtRead(void)
Definition parser.cpp:353
bool nMatch(const char *chKey) const
Definition parser.h:135
int strcmp(const char *s2)
Definition parser.h:177
bool lgEOL(void) const
Definition parser.h:98
NORETURN void NoNumb(const char *chDesc) const
Definition parser.cpp:233
bool m_lgEOF
Definition parser.h:42
t_elementnames elementnames
t_input input
Definition input.cpp:12
t_save save
Definition save.cpp:5
void save_average(long int ipPun)
void parse_save_average(Parser &p, long int ipPun, char *chHeader)
static long int nLine