signal_alloc.c

LastWave

/*..........................................................................*/
/*                                                                          */
/*      L a s t W a v e    P a c k a g e 'signal' 2.1                       */
/*                                                                          */
/*      Copyright (C) 1998-2003 Emmanuel Bacry.                             */
/*      email : [e-mail]                              */
/*                                                                          */
/*..........................................................................*/
/*                                                                          */
/*      This program is a free software, you can redistribute it and/or     */
/*      modify it under the terms of the GNU General Public License as      */
/*      published by the Free Software Foundation; either version 2 of the  */
/*      License, or (at your option) any later version                      */
/*                                                                          */
/*      This program is distributed in the hope that it will be useful,     */
/*      but WITHOUT ANY WARRANTY; without even the implied warranty of      */
/*      MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       */
/*      GNU General Public License for more details.                        */
/*                                                                          */
/*      You should have received a copy of the GNU General Public License   */
/*      along with this program (in a file named COPYRIGHT);                */
/*      if not, write to the Free Software Foundation, Inc.,                */
/*      59 Temple Place, Suite 330, Boston, MA  02111-1307  USA             */
/*                                                                          */
/*..........................................................................*/


/****************************************************************************/
/*                                                                          */
/*  signal_alloc.c   Functions which deal with the dynamical                */
/*                   allocation of memory for SIGNAL's                      */
/*                                                                          */
/****************************************************************************/



#include "lastwave.h"
#include "signals.h"
#include "int_fsilist.h"



/* Default name for a signal */
static char defaultName[] = "";

/* Types of a signal */
char *signalType = "&signal";
char *signaliType = "&signali";


/*
 * Print a signal when it is a result
 */

#define YSigPrintLength 6
#define XYSigPrintLength 3
 
void PrintSignal(SIGNAL signal)
{
  int n;

  if (signal->size == 0) {
    Printf("<size=%d>",signal->size);
    return;
  }
  for (n = 0; n < signal->size; n++) {
    Printf("% 7g ", XSig(signal,n));
    Printf(LWDoubleFormat,signal->Y[n]);
    Printf("\n");
  }
}

/*
 * String conversion
 */

char *ToStrSignal(SIGNAL signal, char flagShort)
{
  static char strShort[50],strShort1[50];
  int n, nAlloc;
  char *str;
  
  if (flagShort) {
    if (signal->type == YSIG) sprintf(strShort,"<size=%d>",signal->size);
    else sprintf(strShort,"<XY,size=%d>",signal->size);
    return(strShort);
  }

  if (signal->size == 0) {
    sprintf(strShort,"<size=%d>",signal->size);
    return(strShort);
  }

  nAlloc = 300;
  str = CharAlloc(nAlloc);
  TempPtr(str);
  if (signal->type == YSIG) {
    if (signal->size < YSigPrintLength) sprintf(str,"<");
    else {
      sprintf(str,"<size=%d",signal->size);
      if (signal->dx == 1 && signal->x0 == 0) strcat(str,";");
      else strcat(str,",");
    }
    if (signal->dx != 1 || signal->x0 != 0) {
      sprintf(strShort,"x0=%g,dx=%g;",signal->x0,signal->dx);
      strcat(str,strShort);
    }
    for (n=0;n<MIN(YSigPrintLength,signal->size);n++) {
      if (strlen(str) > nAlloc-40) {
        nAlloc+=300;
        str = CharAlloc(nAlloc);
        TempPtr(str);
      }
      sprintf(strShort,LWDoubleFormat,signal->Y[n]);
      strcat(str,strShort);
      if (n!=MIN(YSigPrintLength,signal->size)-1) strcat(str,",");
    }
    if (signal->size > YSigPrintLength) strcat(str,",...>");
    else strcat(str,">");
    return(str);
  }
  else {
    if (signal->size < XYSigPrintLength) sprintf(str,"<");
    else sprintf(str,"<size=%d;",signal->size);
    for (n=0;n<MIN(XYSigPrintLength,signal->size);n++) {
      if (strlen(str) > nAlloc-40) {
        nAlloc+=300;
        str = CharAlloc(nAlloc);
        TempPtr(str);
      }
      sprintf(strShort,"(");
      sprintf(strShort1,LWDoubleFormat,signal->X[n]);
      strcat(strShort,strShort1);
      strcat(strShort,"/");
      sprintf(strShort1,LWDoubleFormat,signal->Y[n]);
      strcat(strShort,strShort1);
      strcat(strShort,")");
      strcat(str,strShort);
      if (n!=MIN(XYSigPrintLength,signal->size)-1) strcat(str,",");
    }
    if (signal->size > XYSigPrintLength) strcat(str,",...>");
    else strcat(str,">");
    return(str);
  }
}

/*
 * Print the info of a signal
 */
void PrintInfoSignal(SIGNAL signal)
{
  if (signal->type == XYSIG) Printf("   XY Signal\n");
  else  Printf("   Y Signal\n");

  Printf("   size    :  %d  [%dx %dy]\n",signal->size,signal->sizeMallocX,signal->sizeMallocY);
  
  if (signal->type == YSIG) {
    Printf("   x0      :  %g \n",signal->x0);
    Printf("   dx      :  %g \n",signal->dx);
  }

  Printf("   firstp  :  %d (%g)\n",signal->firstp,XSig(signal,signal->firstp));
  Printf("   lastp   :  %d (%g)\n",signal->lastp,XSig(signal,signal->lastp));
}


/*
 * The different extract options 
 */
static char *extractOptionsSig[] = {"*nolimit","*bperiodic","*bmirror","*bmirror1","*bconst","*b0","*x","*xlin",NULL};
static char *indexextractOptions[] = {"*nolimit", NULL};
enum {
  FSIOptSigNoLimit = FSIOption1,
  FSIOptSigBPer = FSIOption2,
  FSIOptSigBMir = FSIOption3,
  FSIOptSigBMir1 = FSIOption4,
  FSIOptSigBCon = FSIOption5,
  FSIOptSigB0 = FSIOption6,
  FSIOptSigX = FSIOption7,
  FSIOptSigXLin = FSIOption8
};


/* 
 * Set a signal using explicit parameters
 */
char *SetSignalField_(SIGNAL sigLeft,char *field, FSIList *fsiList, LWFLOAT fltRight, VALUE value,char *equal, char *fieldName)
{
  FSI_DECL;
  int i,j,k,size,_j,_iold;
  char *type;
  SIGNAL sigRight,sig;
  LISTV lv;
  LWFLOAT *arrayLeft;
  char flagY;
  RANGE rg;
  
  /* Get the type of the right value  */
  if (value == NULL) type = numType;
  else {
    value = ValueOf(value);
    type = GetTypeValue(value);
  }

  if (type == signaliType && field != NULL && !strcmp("X",field) && fsiList == NULL) {
    sigRight = CastValue(value,SIGNAL);
    if (sigLeft->size != sigRight->size) {
      SetErrorf("Size on both sides should be the same (left size = %d and right size = %d)",sigLeft->size,sigRight->size);
      return(NULL);
    }
    if (sigLeft->type == YSIG) SizeSignal(sigLeft,sigLeft->size,XYSIG);
    for (i=0;i<sigLeft->size;i++) sigLeft->X[i] = sigRight->Y[i];
    return(signaliType);
  }
  
  /* Conversion of range to signal */
  if (type == rangeType && (field == NULL || strcmp("X",field) || fsiList != NULL)) {
    rg = CastValue(value,RANGE);
    sig = TNewSignal();
    SizeSignal(sig,rg->size,YSIG);
    sig->Y[0] = rg->first;
    for (i=1;i<rg->size;i++) sig->Y[i] = sig->Y[i-1]+rg->step;
    value = (VALUE) sig;
    type = signaliType;
  }
  
  /* Case s.X = range */
  if (type == rangeType) {
    rg = CastValue(value,RANGE);
    if (rg->size != sigLeft->size) {
      SetErrorf("Size on both sides should be the same (left size = %d and right size = %d)",sigLeft->size,rg->size);
      return(NULL);
    }
    if (rg->step <= 0) {
      SetErrorf("Step of right handside range should be strictly positive");
      return(NULL);
    }
    rg = CastValue(value,RANGE);
    if (sigLeft->type == XYSIG) sigLeft->type = YSIG;
    sigLeft->dx = rg->step;
    sigLeft->x0 = rg->first;
    return(signaliType);
  }
  
  if (type != numType && type != signaliType && type != signalType && type != listvType) {
    SetErrorf("Expect a signal (and not a '%s') on right handside",type);
    return(NULL);
  }

  if (type == numType) {
    if (value != NULL) fltRight = CastValue(value,NUMVALUE)->f;
    sigRight = NULL;
  }
  else if (type == signaliType || type == signalType) {
    sigRight = CastValue(value,SIGNAL);
  }
  else lv = CastValue(value,LISTV);

  
  /*************************
   *
   * Case value operator is *= += -= /= ^= %=
   *
   *************************/
   
  if (*equal != '=' && *equal != ':') {

    if (type == listvType) {
      SetErrorf("Expect a signal or a number (and not a '%s') on right handside",type);
      return(NULL);
    }

    if (sigLeft->size  == 0) {
      SetErrorf("Expect a non empty signal on the left handside");
      return(NULL);
    }
    
    if (sigRight && (fsiList && sigRight->size != fsiList->nx || !fsiList && sigRight->size != sigLeft->size)) {
      SetErrorf("Signals should have the same size");
      return(NULL);
    }

    if (field != NULL && !strcmp(field,"X")) {
      if (sigLeft->type == YSIG) {
        SetErrorf("Sorry, this functionality (for YSIG) is not implemented yet");
        return(NULL);
      }
      arrayLeft = sigLeft->X;
    }
    else arrayLeft = sigLeft->Y;
    
    switch (*equal) {
    case '+' :  
      if (fsiList == NULL || fsiList->nx == 0) {
        if (sigRight == NULL) for (i=0;i<sigLeft->size;i++) arrayLeft[i]+=fltRight;
        else for (i=0;i<sigLeft->size;i++) arrayLeft[i]+=sigRight->Y[i];
      }
      else {
        if (sigRight == NULL) {FSI_FOR_START(fsiList);arrayLeft[_i]+=fltRight;FSI_FOR_END;}
        else {FSI_FOR_START(fsiList); arrayLeft[_i]+=sigRight->Y[_k];FSI_FOR_END}
      }      
      break; 
    case '-' :  
      if (fsiList == NULL || fsiList->nx == 0) {
        if (sigRight == NULL) for (i=0;i<sigLeft->size;i++) arrayLeft[i]-=fltRight;
        else for (i=0;i<sigLeft->size;i++) arrayLeft[i]-=sigRight->Y[i];
      }
      else {
        if (sigRight == NULL) {FSI_FOR_START(fsiList);arrayLeft[_i]-=fltRight;FSI_FOR_END;}
        else {FSI_FOR_START(fsiList); arrayLeft[_i]-=sigRight->Y[_k];FSI_FOR_END}
      }      
      break; 
    case '*' :  
      if (fsiList == NULL || fsiList->nx == 0) {
        if (sigRight == NULL) for (i=0;i<sigLeft->size;i++) arrayLeft[i]*=fltRight;
        else for (i=0;i<sigLeft->size;i++) arrayLeft[i]*=sigRight->Y[i];
      }
      else {
        if (sigRight == NULL) {FSI_FOR_START(fsiList);arrayLeft[_i]*=fltRight;FSI_FOR_END;}
        else {FSI_FOR_START(fsiList); arrayLeft[_i]*=sigRight->Y[_k];FSI_FOR_END}
      }      
      break; 
    case '^' :  
      if (fsiList == NULL || fsiList->nx == 0) {
        if (sigRight == NULL) for (i=0;i<sigLeft->size;i++) arrayLeft[i]=pow(fabs(arrayLeft[i]),fltRight);
        else for (i=0;i<sigLeft->size;i++) arrayLeft[i]=pow(fabs(arrayLeft[i]),sigRight->Y[i]);
      }
      else {
        if (sigRight == NULL) {FSI_FOR_START(fsiList);arrayLeft[_i]=pow(fabs(arrayLeft[_i]),fltRight);FSI_FOR_END;}
        else {FSI_FOR_START(fsiList); arrayLeft[_i]=pow(fabs(arrayLeft[_i]),sigRight->Y[_k]);FSI_FOR_END}
      }      
      break; 
    case '%' :  
      if (fsiList == NULL || fsiList->nx == 0) {
        if (sigRight == NULL) for (i=0;i<sigLeft->size;i++) arrayLeft[i]=((int) arrayLeft[i]) % ((int) fltRight);
        else for (i=0;i<sigLeft->size;i++) arrayLeft[i]=((int) arrayLeft[i]) % ((int) sigRight->Y[i]);
      }
      else {
        if (sigRight == NULL) {FSI_FOR_START(fsiList);arrayLeft[_i]=((int) arrayLeft[_i]) % ((int) fltRight);FSI_FOR_END;}
        else {FSI_FOR_START(fsiList); arrayLeft[_i]=((int) arrayLeft[_i]) % ((int) sigRight->Y[_k]);FSI_FOR_END}
      }      
      break; 
    case '/' :  
      if (sigRight == NULL && fltRight == 0) {      
        SetErrorf("Division by 0");
        return(NULL);
      }
      if (fsiList == NULL || fsiList->nx == 0) {
        if (sigRight == NULL) for (i=0;i<sigLeft->size;i++) arrayLeft[i]/=fltRight;
        else for (i=0;i<sigLeft->size;i++) {
          if (sigRight->Y[i] == 0) {     
            SetErrorf("Division by 0");
            return(NULL);
          }
          arrayLeft[i]/=sigRight->Y[i];
        }
      }
      else {
        if (sigRight == NULL) {FSI_FOR_START(fsiList);arrayLeft[_i]/=fltRight;FSI_FOR_END;}
        else {
          FSI_FOR_START(fsiList); 
          if (sigRight->Y[_k] == 0) {     
            SetErrorf("Division by 0");
            return(NULL);
          }
          arrayLeft[_i]/=sigRight->Y[_k];
          FSI_FOR_END;
        }
      }    
      break; 
    default :   break;
  	}
  	return(signaliType);
  }


  /*************************
   *
   * Case s[rang1...rangeN] = sig or LWFLOAT
   *
   *************************/

  if (*equal == '=' && type != listvType) {
    /* Case of s=sigRight */
    if (fsiList == NULL) {
      if (sigRight) CopySig(sigRight,sigLeft);
      else {
        SizeSignal(sigLeft,1,YSIG);
        sigLeft->Y[0] = fltRight;
      }
      return(signaliType);
    }      
      
    /* Case of s[]=sigRight */
    if (fsiList->nx==0) {
      if (sigRight) CopySig(sigRight,sigLeft);
      else {
        SizeSignal(sigLeft,1,YSIG);
        sigLeft->Y[0] = fltRight;
      }
      return(signaliType);
    }      
    
    if (sigLeft->size  == 0) {
      SetErrorf("Expect a non empty signal on the left handside");
      return(NULL);
    }
   
    if (field != NULL && !strcmp(field,"X")) {
      if (sigLeft->type == YSIG) {
        SetErrorf("Sorry, this functionality (for YSIG) is not implemented yet");
        return(NULL);
      }
      arrayLeft = sigLeft->X;
    }
    else arrayLeft = sigLeft->Y;
  
    /* Case s[1]=2 */
    if (type == numType) {
      if (fsiList->nx != 1) {
        SetErrorf("Size of both handsides should match (left size = %d, right size =1)",fsiList->nx);
        return(NULL);
      }
      arrayLeft[(int) FSIArray((&(fsiList->fsi[0])),0)]=fltRight;
      return(signaliType);
    }
    

    /* Case of s[...]=sigRight */
    if (sigRight->size != fsiList->nx) {
      SetErrorf("Size of both handsides should match (left size = %d, right size = %d)",fsiList->nx,sigRight->size);
      return(NULL);
    }
    FSI_FOR_START(fsiList); 
    arrayLeft[_i] = sigRight->Y[_k];
    FSI_FOR_END;
    return(signaliType);
  }

  /*************************
   *
   * Case of s[range1,...,rangeN] = {sig1,...,sigN} or s[range1,...,rangeN] := {sig}
   *
   *************************/

  if (type == listvType) {
    
    if (fsiList == NULL) {
      SetErrorf("Cannot set a signal with a listv");
      return(NULL);
    }
      
    if (*equal == '=' && lv->length != fsiList->size) {
      SetErrorf("Length of extraction list should be the same as length of &listv");
      return(NULL);
    }

    if (*equal == ':' && lv->length != 1) {
      SetErrorf("Length of right handside &listv is expected to be 1");
      return(NULL);
    }
    
    if (field != NULL && !strcmp(field,"X")) {
      SetErrorf("Sorry, this functionality is not implemented yet");
      return(NULL);
    }
    
    if (field != NULL && !strcmp(field,"Y")) flagY = YES;
    else flagY = NO;
       
    /* We first need to compute the new size */
    size = 0;
    for (j=0;j<lv->length;j++) {
      if (lv->values[j] != NULL && GetTypeValue(lv->values[j]) != signaliType && GetTypeValue(lv->values[j]) != signalType ) {
        SetErrorf("The &listv should contain only numbers or signals");
        return(NULL);
      }
      size += (lv->values[j] == NULL ? 1 : CastValue(lv->values[j],SIGNAL)->size);
    }
    if (*equal == ':') size *= fsiList->size;
    
    /* Then we go */
    size = sigLeft->size + size -fsiList->nx;
    sig = NewSignal();
    SizeSignal(sig,size,sigLeft->type);
  
    /* Then we perform the set */
    j = 0; /* Index in the sigLeft string */
    k = 0; /* Index in the sig string */
    _iold = -1;
    FSI_FOR_START1(fsiList);
    if (*equal == ':') GetListvNth(lv,0,(VALUE *) &sigRight,&fltRight);
    else GetListvNth(lv,_n,(VALUE *) &sigRight,&fltRight);
    switch(fsi->type) {
    case FSIFloat : 
      _j = (int) fsi->val.f;
      if (_j<=_iold) {
        SetErrorf("Sorry, this operation with non sorted indexes is not implemented");
        DeleteSignal(sig);
        return(NULL);
      }
      _iold = _j;
      memcpy(sig->Y+k,sigLeft->Y+j,(_j-j)*sizeof(LWFLOAT));
      if (!flagY && sigLeft->type == XYSIG) memcpy(sig->X+k,sigLeft->X+j,(_j-j)*sizeof(LWFLOAT));
      k+=_j-j;
      if (sigRight) {
        if (sigRight->size != 0) {
          memcpy(sig->Y+k,sigRight->Y,sizeof(LWFLOAT)*sigRight->size);
          if (!flagY && sig->type == XYSIG) {
            if (sigRight->type == XYSIG) memcpy(sig->X+k,sigRight->X,sizeof(LWFLOAT)*sigRight->size);
            else for (i=0;i<sigRight->size;i++) sig->X[k+i]= sigLeft->X[(_j-1<0 ? 0 : _j-1)];
          }
          k+=sigRight->size;
        }
      }
      else {sig->Y[k]=fltRight;k++;}
      j=_j+1;
      break;
    case FSIRange : 
      _j = (int) fsi->val.r->first;
      if (_j<=_iold) {
        SetErrorf("Sorry, this operation with non sorted indexes is not implemented");
        DeleteSignal(sig);
        return(NULL);
      }