I'm trying to set up an integration of spectral emission over a sphere. The problem is that the integral has to be done over longitude and latitude "theta" and "phi", but I need a spectrum to come out as a function of wavelength, "lambda". Wavelength is NOT a variable of integration -- rather, an array of wavelengths goes in and an array of fluxes as a function of wavelength comes out.
Right now I have this implemented as a for loop over wavelength inside which are the two nested "qpint1d" function calls. The first "qpint1d" function also has a for loop in it over "theta" to make the inner "qpint1d" call work with "theta" as a scalar. It works but as you can imagine, with nested for loops this is not very efficient.
I've tried to pass the wavelength through to the inner integral as a vector within the "private" structure so that instead of a vector over "phi", the inner function would return a matrix of "phi" vs. wavelength, thus eliminating the outermost for loop (over wavelength). No luck, though -- I'm just getting the error message from QPINT1D_QKEVAL about how the integrand function must return a vector of values. Is there a better way to do this?
On Tuesday, November 13, 2012 1:28:16 PM UTC-5, lucylim wrote:
> Hi all,
> I'm trying to set up an integration of spectral emission over a sphere. The problem is that the integral has to be done over longitude and latitude "theta" and "phi", but I need a spectrum to come out as a function of wavelength, "lambda". Wavelength is NOT a variable of integration -- rather, an array of wavelengths goes in and an array of fluxes as a function of wavelength comes out.
> Right now I have this implemented as a for loop over wavelength inside which are the two nested "qpint1d" function calls. The first "qpint1d" function also has a for loop in it over "theta" to make the inner "qpint1d" call work with "theta" as a scalar. It works but as you can imagine, with nested for loops this is not very efficient.
> I've tried to pass the wavelength through to the inner integral as a vector within the "private" structure so that instead of a vector over "phi", the inner function would return a matrix of "phi" vs. wavelength, thus eliminating the outermost for loop (over wavelength). No luck, though -- I'm just getting the error message from QPINT1D_QKEVAL about how the integrand function must return a vector of values. Is there a better way to do this?
> Many thanks,
> Lucy Lim
> NASA/GSFC
You should just come down the hall and ask :-)
The only vector input to your integrand function should be the variable of integration. You are right the the wavelength part will need to be done as a FOR loop. Sorry.
And yes, you will need to do an inner loop over THETA variable. The "1D" of QPINT1D is there for a reason. Sorry. Since QPINT1D does a fair amount of work per call, I don't think it would be any faster to try to "vectorize" it.
On Tuesday, November 13, 2012 4:42:39 PM UTC-5, Craig Markwardt wrote:
> On Tuesday, November 13, 2012 1:28:16 PM UTC-5, lucylim wrote:
> > Hi all,
> > I'm trying to set up an integration of spectral emission over a sphere. The problem is that the integral has to be done over longitude and latitude "theta" and "phi", but I need a spectrum to come out as a function of wavelength, "lambda". Wavelength is NOT a variable of integration -- rather, an array of wavelengths goes in and an array of fluxes as a function of wavelength comes out.
> > Right now I have this implemented as a for loop over wavelength inside which are the two nested "qpint1d" function calls. The first "qpint1d" function also has a for loop in it over "theta" to make the inner "qpint1d" call work with "theta" as a scalar. It works but as you can imagine, with nested for loops this is not very efficient.
> > I've tried to pass the wavelength through to the inner integral as a vector within the "private" structure so that instead of a vector over "phi", the inner function would return a matrix of "phi" vs. wavelength, thus eliminating the outermost for loop (over wavelength). No luck, though -- I'm just getting the error message from QPINT1D_QKEVAL about how the integrand function must return a vector of values. Is there a better way to do this?
> > Many thanks,
> > Lucy Lim
> > NASA/GSFC
> You should just come down the hall and ask :-)
> The only vector input to your integrand function should be the variable of integration. You are right the the wavelength part will need to be done as a FOR loop. Sorry.
> And yes, you will need to do an inner loop over THETA variable. The "1D" of QPINT1D is there for a reason. Sorry. Since QPINT1D does a fair amount of work per call, I don't think it would be any faster to try to "vectorize" it.
> Craig
OK -- nested FOR loops it is, then. (Just don't tell JD.) Thanks! Also, it's immensely helpful to have the quadpack available in IDL, so thanks for writing that too.
