| GAUSS_FUNCT problem [message #38230] |
Fri, 27 February 2004 09:11  |
Michael Wallace
Messages: 409
Registered: December 2003
|
Senior Member |
|
|
I have discovered something interesting in the RSI-provided
gaussfit.pro. I hesitate to call this a "bug," but I don't know what
else to call it. Hopefully, it's just something stupid I'm doing. I
have copied the pertinent documentation and code below.
; CALLING SEQUENCE:
; FUNCT,X,A,F,PDER
; INPUTS:
; X = VALUES OF INDEPENDENT VARIABLE.
; A = PARAMETERS OF EQUATION DESCRIBED BELOW.
; OUTPUTS:
; F = VALUE OF FUNCTION AT EACH X(I).
;
; OPTIONAL OUTPUT PARAMETERS:
; PDER = (N_ELEMENTS(X),6) ARRAY CONTAINING THE
; PARTIAL DERIVATIVES. P(I,J) = DERIVATIVE
; AT ITH POINT W/RESPECT TO JTH PARAMETER.
; PROCEDURE:
; F = A(0)*EXP(-Z^2/2) + A(3) + A(4)*X + A(5)*X^2
; Z = (X-A(1))/A(2)
; Elements beyond A(2) are optional.
PRO GAUSS_FUNCT,X,A,F,PDER
COMPILE_OPT idl2, hidden
ON_ERROR,2 ;Return to caller if an error occurs
n = n_elements(a)
if a[2] ne 0.0 then begin
Z = (X-A[1])/A[2] ;GET Z
EZ = EXP(-Z^2/2.) ;GAUSSIAN PART
endif else begin
z = 100.
ez = 0.0
endelse
case n of
3: F = A[0]*EZ
4: F = A[0]*EZ + A[3]
5: F = A[0]*EZ + A[3] + A[4]*X
6: F = A[0]*EZ + A[3] + A[4]*X + A[5]*X^2 ;FUNCTIONS.
ENDCASE
The variable X is an array of values. In the case where A[2] is not
equal to 0, Z and EZ are created as arrays. However, when A[2] is equal
to 0, Z and EZ are created as scalars. Then the value of F is computed
by using EZ. If EZ is an array, F is an array. If EZ is scalar, F is
scalar. F should be the value of the function of each X[I]. F should
be an array the same size as X, but when A[2] is equal to 0, F is just a
single scalar value. The problem is that by the contract of procedure
states that F will be an array. And GAUSSFIT itself expects an array,
and this causes problems!!
Of course, the "fix" to this is to make Z an array of N elements with
each element set to 100 and EZ an array of N elements with each element
set to 0 in the case where A[2] is equal to 0. This ensures that F is
always an array.
So, is this a real error that needs to get fixed or am I imagining things?
-Mike
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
