closes #1843 --- interpol supports complex values. (#1844)

gnudatalanguage · Jun 12, 2024 · bb6d739 · bb6d739
1 parent d2978d1
commit bb6d739
Show file tree

Hide file tree

Showing 2 changed files with 24 additions and 4 deletions.
diff --git a/src/interpol.cpp b/src/interpol.cpp
@@ -374,6 +374,8 @@ namespace lib {
     //dimensions
     BaseGDL* p0 = e->GetParDefined(0);
     DType type=p0->Type();
+	if (type==GDL_STRING) e->Throw("Internal Error, string value not allowed when calling GDL_INTERPOL, please report.");
+	if (ComplexType(type)) e->Throw("Internal Error, complex value not allowed when calling GDL_INTERPOL, please report.");
     SizeT nV=p0->N_Elements();
     SizeT orig_nV=nV;
     if (nV <nmin) e->Throw("V has too few elements for this kind of interpolation.");
@@ -395,7 +397,7 @@ namespace lib {
     //type of output
     bool isDouble = false;
 
-    if (nParam==2) {       // we make X, only V has to be chacked for Nan
+    if (nParam==2) {       // we make X, only V has to be checked for Nan
       BaseGDL* p1 = e->GetParDefined(1);
       if (p1->N_Elements() >1) e->Throw("N must be one positive integer");
       DLongGDL* n1gdl=e->GetParAs<DLongGDL>(1);

diff --git a/src/pro/utilities/interpol.pro b/src/pro/utilities/interpol.pro
@@ -7,12 +7,30 @@
 function INTERPOL, p0, p1, p2, _EXTRA=extra
 ;
 COMPILE_OPT idl2, HIDDEN
-ON_ERROR, 2
+ON_ERROR, 1
 ;
+; test for valid params
+; P0 (aka V) must not be a string
+if isa(P0,/STRING) then Message,"Operation illegal with strings."
+; sort arrays if necessary
 if n_params() eq 3 then begin
    s=sort(p1)
    p1=p1[s]
    p0=p0[s]
-   return,GDL_INTERPOL(p0, p1, p2, _EXTRA=extra)
-endif else return,GDL_INTERPOL( p0, p1, _EXTRA=extra)
+endif
+; treat complex types as two real types, below
+if isa(P0,/COMPLEX) then begin
+    if n_params() eq 3 then begin
+       real=GDL_INTERPOL(real_part(p0), p1, p2, _EXTRA=extra)
+       imag=GDL_INTERPOL(imaginary(p0), p1, p2, _EXTRA=extra)
+       return,complex(real,imag)
+    endif else begin
+        real=GDL_INTERPOL(real_part(p0), p1, _EXTRA=extra)
+        imag=GDL_INTERPOL(imaginary(p0), p1, _EXTRA=extra)
+        return,complex(real,imag)
+    endelse
+endif
+
+; all others
+if n_params() eq 3 then return,GDL_INTERPOL(p0, p1, p2, _EXTRA=extra) else return,GDL_INTERPOL( p0, p1, _EXTRA=extra)
 end