cloudy  trunk
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
cont_ipoint.cpp
Go to the documentation of this file.
1 /* This file is part of Cloudy and is copyright (C)1978-2017 by Gary J. Ferland and
2  * others. For conditions of distribution and use see copyright notice in license.txt */
3 /*ipoint returns pointer to any energy within energy mesh */
4 /*ipFineCont returns array index within fine energy mesh */
5 /*ipContEnergy generate pointer to energy within continuum array
6  * continuum energy in Rydbergs */
7 /*ipLineEnergy generate pointer to line energy within energy mesh
8  * line energy in Rydbergs */
9 #include "cddefines.h"
10 #include "ipoint.h"
11 #include "prt.h"
12 #include "rfield.h"
13 
14 /*ipoint returns pointer to any energy within energy mesh */
15 long ipoint(double energy_ryd)
16 {
17  if( energy_ryd < rfield.emm() || energy_ryd > rfield.egamry() )
18  {
19  fprintf( ioQQQ, " ipoint:\n" );
20  fprintf( ioQQQ, " The energy_ryd array is not defined at nu=%11.3e. The bounds are%11.3e%11.3e\n",
21  energy_ryd, rfield.emm(), rfield.egamry() );
22  fprintf( ioQQQ, " ipoint is aborting to get trace, to find how this happened\n" );
23  ShowMe();
25  }
26 
27  return rfield.ipointF(energy_ryd);
28 }
29 
30 /*ipContEnergy generate pointer to energy within continuum array */
32  /* continuum energy in Rydbergs */
33  double energy,
34  /* 4 char label for continuum, like those returned by chLineLbl */
35  const char *chLabel)
36 {
37  long int ipConSafe_v;
38 
39  DEBUG_ENTRY( "ipContEnergy()" );
40 
41  ipConSafe_v = ipoint(energy);
42 
43  /* write label in this cell if not anything there yet */
44  if( strcmp(rfield.chContLabel[ipConSafe_v-1].c_str()," ") == 0 )
45  {
46  rfield.chContLabel[ipConSafe_v-1] = chLabel ;
47  }
48 
49  /* this is a quick way to find all continua that occur within a given cell,
50  * recall the off-by-one aspect of C */
51  {
52  enum {DEBUG_LOC=false};
53  if( DEBUG_LOC )
54  {
55  /* recall the off-by-one aspect of C - the number below is
56  * on the physics scale because this returns the index
57  * on that scale, so the first is 1 for [0] */
58  if( ipConSafe_v == 23 )
59  fprintf(ioQQQ,"%s\n", chLabel );
60  }
61  }
62  return ipConSafe_v;
63 }
64 
65 /*ipLineEnergy generate pointer to line energy within energy mesh
66  * line energy in Rydbergs -
67  * return value is array index on the physics or Fortran scale, counting from 1 */
68 long ipLineEnergy(double energy,
69  /* 4 char label for line, like those returned by chLineLbl */
70  const char *chLabel ,
71  /* will make sure energy is < this array index if greater than 0, does nothing if <= 0*/
72  long ipIonEnergy )
73 {
74  long int ipLine_ret;
75 
76  DEBUG_ENTRY( "ipLineEnergy()" );
77 
78  ipLine_ret = ipoint(energy);
79  ASSERT( ipLine_ret );
80  /* make sure pointer is below next higher continuum if positive */
81  if( ipIonEnergy > 0 )
82  {
83  ipLine_ret = MIN2( ipLine_ret , ipIonEnergy-1 );
84  }
85 
86  ASSERT( ipLine_ret > 0 );
87  /* stuff in a label if none there,
88  * note that this is offset one below actual number to be on C scale of arrays */
89  /* >>chng 06 nov 23, faster to use line_count index rather than checking 5 chars
90  * first call will have zero so false */
91  /*if( strcmp(rfield.chLineLabel[ipLine_ret-1]," ")==0 )*/
92  if( !rfield.line_count[ipLine_ret-1] )
93  {
94  rfield.chLineLabel[ipLine_ret-1] = chLabel;
95  }
96  /* this keeps track of the number of lines within this cell */
97  ++rfield.line_count[ipLine_ret-1];
98 
99  /* this is a quick way to find all lines that occur within a given cell,
100  * recall the off-by-one aspect of C */
101  {
102  enum {DEBUG_LOC=false};
103  if( DEBUG_LOC )
104  {
105  /* recall the off-by-one aspect of C - the numbers is one more
106  * than the index on the C scale */
107  if( ipLine_ret == 23 )
108  fprintf(ioQQQ,"%s\n", chLabel );
109  }
110  }
111 
112  /* this implements the print continuum indices command */
113  if( prt.lgPrtContIndices )
114  {
115  /* print header if first time */
116  static bool lgFirst = true;
117  if( lgFirst )
118  {
119  /* print header and set flag so do not do again */
120  fprintf(ioQQQ , "\n\noutput from print continuum indices command follows.\n");
121  fprintf(ioQQQ , "cont ind (F scale)\tenergy(ryd)\tlabel\n");
122  lgFirst = false;
123  }
124  if( energy >= prt.lgPrtContIndices_lo_E && energy <= prt.lgPrtContIndices_hi_E )
125  {
126  /* use varying formats depending on order of magnitude of energy
127  * NB - the ipLine_ret is the index on the physical or Fortran scale,
128  * and is 1 greater than the array element used in the code, which is
129  * on the C scale */
130  if( energy < 1. )
131  {
132  fprintf(ioQQQ , "%li\t%.3e\t%s\n" , ipLine_ret , energy , chLabel);
133  }
134  else if( energy < 10. )
135  {
136  fprintf(ioQQQ , "%li\t%.3f\t%s\n" , ipLine_ret , energy , chLabel);
137  }
138  else if( energy < 100. )
139  {
140  fprintf(ioQQQ , "%li\t%.2f\t%s\n" , ipLine_ret , energy , chLabel);
141  }
142  else
143  {
144  fprintf(ioQQQ , "%li\t%.1f\t%s\n" , ipLine_ret , energy , chLabel);
145  }
146  }
147  }
148 
149  if( prt.lgPrnLineCell )
150  {
151  /* print line cell option - number on command line is cell on Physics scale */
152  if( prt.nPrnLineCell == ipLine_ret )
153  {
154  static bool lgMustPrintHeader = true;
155  if( lgMustPrintHeader )
156  fprintf(ioQQQ, "Lines within cell %li (physics scale) \nLabel\tEnergy(Ryd)\n",prt.nPrnLineCell );
157  lgMustPrintHeader = false;
158  fprintf(ioQQQ,"%s\t%.3e\n" , chLabel , energy );
159  }
160  }
161  return ipLine_ret;
162 }
163 
164 /*ipFineCont returns array index within fine energy mesh */
166  /* energy in Ryd */
167  double energy_ryd )
168 {
169  long int ipoint_v;
170 
171  DEBUG_ENTRY( "ipFineCont()" );
172 
173  if( energy_ryd < rfield.fine_ener_lo || energy_ryd > rfield.fine_ener_hi )
174  {
175  return -1;
176  }
177 
178  /* this is on the Fortran scale of array index counting, so is one above the
179  * c scale. later the -1 will appear on any index references
180  *
181  * 07 Jun 22 testing done to confirm that energy grid is correct: did
182  * same sim with standard fine continuum resolution, and another with 200x
183  * finer resolution, and confirmed that these line up correctly. */
184  ipoint_v = (long int)(log10(energy_ryd*(1.-rfield.fine_resol/2.) /
185  rfield.fine_ener_lo)/log10(1.+rfield.fine_resol));
186 
187  ASSERT( ipoint_v >= 0 && ipoint_v< rfield.nfine );
188  return ipoint_v;
189 }
double emm() const
Definition: mesh.h:84
long int * line_count
Definition: rfield.h:61
realnum lgPrtContIndices_hi_E
Definition: prt.h:203
vector< string > chContLabel
Definition: rfield.h:215
long int nPrnLineCell
Definition: prt.h:256
long ipFineCont(double energy_ryd)
FILE * ioQQQ
Definition: cddefines.cpp:7
#define MIN2(a, b)
Definition: cddefines.h:807
double fine_resol
Definition: rfield.h:387
long ipoint(double energy_ryd)
Definition: cont_ipoint.cpp:15
bool lgPrnLineCell
Definition: prt.h:253
long ipContEnergy(double energy, const char *chLabel)
Definition: cont_ipoint.cpp:31
double energy(const genericState &gs)
bool lgPrtContIndices
Definition: prt.h:199
long ipLineEnergy(double energy, const char *chLabel, long ipIonEnergy)
Definition: cont_ipoint.cpp:68
t_rfield rfield
Definition: rfield.cpp:9
#define EXIT_FAILURE
Definition: cddefines.h:168
#define cdEXIT(FAIL)
Definition: cddefines.h:484
size_t ipointF(double anu) const
Definition: mesh.h:159
static bool lgMustPrintHeader
Definition: save_line.cpp:264
vector< string > chLineLabel
Definition: rfield.h:212
t_prt prt
Definition: prt.cpp:14
realnum fine_ener_lo
Definition: rfield.h:383
#define ASSERT(exp)
Definition: cddefines.h:617
realnum fine_ener_hi
Definition: rfield.h:383
long nfine
Definition: rfield.h:385
#define DEBUG_ENTRY(funcname)
Definition: cddefines.h:729
double egamry() const
Definition: mesh.h:88
int fprintf(const Output &stream, const char *format,...)
Definition: service.cpp:1217
realnum lgPrtContIndices_lo_E
Definition: prt.h:203
void ShowMe(void)
Definition: service.cpp:205