ecat_roi.c

Go to the documentation of this file.

 
/*****************************************************************************/
#include "libtpcroi.h"
/*****************************************************************************/
 
/*****************************************************************************/
char roierrmsg[128];
/*****************************************************************************/
#define SUPPORT_LEGACY
/*****************************************************************************/
 
/*****************************************************************************/
void roi_empty(
  ROI_list *rl
) {
  int r;
  RoiList *next;
#ifdef SUPPORT_LEGACY 
  if(rl->roi) {
    for(r=0; r<rl->nr; r++) {
      if(rl->roi[r].point_nr>0) {
        free(rl->roi[r].x);
        free(rl->roi[r].y);
        rl->roi[r].point_nr=0;
      }
    }
  }
  if(rl->_allocNr>0) free(rl->roi);
  rl->_allocNr=0;
#endif
  rl->nr=0;
  while(rl->rois) {
    next=rl->rois->next;
    free(rl->rois);
    rl->rois=next;
  }
}
/*****************************************************************************/
 
/*****************************************************************************/
void roi_init(
  ROI_list *rl
) {
  rl->roi=(ROI*)NULL;
  rl->nr=0;
#ifdef SUPPORT_LEGACY
  rl->_allocNr=0;
  rl->rois=NULL;
#endif
}
/*****************************************************************************/
 
/*****************************************************************************/
void roi_delete_n(
  ROI_list *rl, 
  int ind
) {
  RoiList *rois=rl->rois;
  int pos=0;
  while(pos<ind) {
    rois=rois->next;
    if(!rois) return;
    pos++;
  }
  if(rois->prev) rois->prev->next=rois->next; else rl->rois=rois->next;
  if(rois->next) rois->next->prev=rois->prev;
  if(rois->roi->x) free(rois->roi->x);
  if(rois->roi->y) free(rois->roi->y);
  free(rois->roi);
  free(rois);
  rl->nr--;
}
/*****************************************************************************/
 
/*****************************************************************************/
void roi_delete(
  ROI_list *rl, 
  ROI *roi
) {
  RoiList *rois=rl->rois;
  while(rois->roi!=roi) {
    rois=rois->next;
    if(!rois) return;
  }
  if(rois->prev) rois->prev->next=rois->next; else rl->rois=rois->next;
  if(rois->next) rois->next->prev=rois->prev;
  if(rois->roi->x) free(rois->roi->x);
  if(rois->roi->y) free(rois->roi->y);
  free(rois->roi);
  free(rois);
  rl->nr--;
}
/*****************************************************************************/
 
/*****************************************************************************/
int roi_read(
  const char *fname, 
  ROI_list *rl
) {
  RoiList *newroientry;
  RoiList *roilist=rl->rois;
  ROI *newroi;
  int verbose=0;
 
  if(verbose>0) printf("roi_read(%s, rl)\n", fname);
  if(fname==NULL || rl==NULL) return(1);
 
  /* read the lines of the ROI file */
  STR_TOKEN_LIST lines;
  str_token_list_init(&lines);
  if(textfileReadLines(fname, &lines)!=0) {
    strcpy(roierrmsg, "cannot read file");
    str_token_list_empty(&lines); return(1);
  }
  if(verbose>3)
    for(int i=0; i<lines.token_nr; i++)
      printf("%d: %s\n", 1+i, lines.tok[i]);
 
  /* Prepare for the data */
  if(roilist) while(roilist->next) roilist=roilist->next;
  int n=0; // roi number
 
  int li=0;
  while(li<lines.token_nr) {
 
    /* Allocate space for the new ROI and initialize it*/
    newroi=malloc(sizeof(ROI));
    memset(newroi, 0, sizeof(ROI));
    newroi->userdata=NULL;
 
    /* Find first line that starts with '*' */
    for(; li<lines.token_nr; li++) if(lines.tok[li][0]=='*') break;
    if(li==lines.token_nr) break;
 
    /* Read line contents until roiname */
    int c=sscanf(lines.tok[li], "*%s %f %f %d %d %d %d %d %d %d %d %d ",
                 newroi->imgfile, &newroi->zoom, &newroi->recon_zoom,
                 &newroi->matnum, &newroi->type, &newroi->status, 
                 &newroi->pos_x, &newroi->pos_y, &newroi->w, &newroi->h, &newroi->t, 
                 &newroi->roi);
    if(c!=12) {
      strcpy(roierrmsg, "invalid format");
      str_token_list_empty(&lines); return(2);
    }
 
    /* read roi name and the rest of the line */
    char buf[256]="", buf2[256];
    int ti=13;
    while(strTokenNCpy(lines.tok[li], " \t", ti, buf2, 256)>0) {
      if(strlen(buf)>0) strcat(buf, " ");
      strcat(buf, buf2);
      ti++;
    }
    if(strlen(buf)<3) {
      strcpy(roierrmsg, "invalid roi name");
      free(newroi->x); free(newroi->y); free(newroi);
      str_token_list_empty(&lines); return(3);
    }
    char *cptr=strstr(buf, "///"); if(cptr!=NULL) {*cptr=(char)0; cptr+=3;}
    c=strncpyCleanSpaces(newroi->roiname, buf, 128);
    if(c==0) sprintf(newroi->roiname, "ROI%d", 1+n);
    strReplaceChar(newroi->roiname, ' ', '_');
    /* Read point number, skipping useless value before it */
    c=sscanf(cptr, "%*d %d \n", &newroi->point_nr);
    if(c!=1) {
      strcpy(roierrmsg, "invalid roi point number");
      free(newroi->x); free(newroi->y); free(newroi);
      str_token_list_empty(&lines); return(4);
    }
 
    /* Trace ROI should contain another line with points */
    if(newroi->type==ROI_TRACE && newroi->point_nr>0) {
      li++; //printf("  line='%s'\n", lines.tok[li]);
      if(strTokenNr(lines.tok[li], " \t") < 2*newroi->point_nr) {
        strcpy(roierrmsg, "cannot read ROI border definition");
        free(newroi->x); free(newroi->y); free(newroi);
        str_token_list_empty(&lines); return(7);
      }
      newroi->x=malloc(newroi->point_nr*sizeof(int));
      newroi->y=malloc(newroi->point_nr*sizeof(int));
      char xs[16], ys[16];
      int i;
      for(i=0; i<newroi->point_nr; i++) {
        if(strTokenNCpy(lines.tok[li], " \t", 2*i+1,   xs, 16)<1) break;
        if(strTokenNCpy(lines.tok[li], " \t", 2*i+2, ys, 16)<1) break;
        newroi->x[i]=atoi(xs);
        newroi->y[i]=atoi(ys);
      } //printf("i=%d / %d\n", i, newroi->point_nr);
      if(i < newroi->point_nr) {
        strcpy(roierrmsg, "cannot read ROI border definition");
        free(newroi->x); free(newroi->y); free(newroi);
        str_token_list_empty(&lines); return(9);
      }
    }
    /* Append ROI to list */
    newroientry=malloc(sizeof(RoiList));
    newroientry->roi=newroi;
    newroientry->prev=roilist;
    newroientry->next=NULL;
    if(roilist) {
      roilist->next=newroientry;
      roilist=newroientry;
    } else {
      roilist=newroientry;
      rl->rois=roilist;
    }
    rl->nr++;
    n++; li++;
  }
  str_token_list_empty(&lines);
 
  if(n==0) {
    strcpy(roierrmsg, "error in ROI file");
    return(2);
  }
  return(0);
}
/*****************************************************************************/
 
/*****************************************************************************/
int roi_compute_rect(
  ROI *roi
) {
  roi->point_nr=4;
  if((roi->x=malloc(roi->point_nr*sizeof(int)))==NULL)
    return(-1);
  if((roi->y=malloc(roi->point_nr*sizeof(int)))==NULL) {
    free(roi->x); return(-1);}
  roi->x[0]=0; roi->y[0]=0;
  roi->x[1]=roi->w; roi->y[1]=0;
  roi->x[2]=roi->w; roi->y[2]=roi->h;
  roi->x[3]=0; roi->y[3]=roi->h;
  return(0);
}
 
#ifdef SUPPORT_LEGACY
int roiComputeRect(
  ROI *roi
) {
  roi->point_nr=5;
  if((roi->x=malloc(roi->point_nr*sizeof(int)))==NULL)
    return(-1);
  if((roi->y=malloc(roi->point_nr*sizeof(int)))==NULL) {
    free((char*)roi->x); return(-1);}
  roi->x[0]=0; roi->y[0]=0;
  roi->x[1]=roi->w; roi->y[1]=0;
  roi->x[2]=roi->w; roi->y[2]=roi->h;
  roi->x[3]=0; roi->y[3]=roi->h;
  roi->x[4]=0; roi->y[4]=0;
  return(0);
}
#endif
/*****************************************************************************/
 
/*****************************************************************************/
int roi_compute_circle(
  ROI *roi
) {
  int points,r;
  r=roi->w/2;
  roi->pos_x+=r;
  roi->pos_y+=roi->w;
  roi->point_nr=36; // TODO: clever algorithm to decide how many points we need
  roi->x=malloc(sizeof(int)*roi->point_nr);
  roi->y=malloc(sizeof(int)*roi->point_nr);
  for(points=0;points<roi->point_nr;points++) {
    float rad=(M_PI/180.0)*(360.0*(points/(float)roi->point_nr));
    roi->x[points]=sin(rad)*r;
    roi->y[points]=cos(rad)*r-r;
  }
  return 0;
}
/*****************************************************************************/
 
/*****************************************************************************/
int roi_compute_ellipse(
  ROI *roi
) {
  int points;
  roi->pos_x+=roi->w/2;
  roi->pos_y+=roi->h;
  roi->point_nr=40; // TODO: clever algorithm to decide how many points we need
  roi->x=malloc(sizeof(int)*roi->point_nr);
  roi->y=malloc(sizeof(int)*roi->point_nr);
  for(points=0;points<roi->point_nr;points++) {
    float rad=(M_PI/180.0)*(360.0*(points/(float)roi->point_nr));
    roi->x[points]=sin(rad)*(roi->w/2);
    roi->y[points]=cos(rad)*(roi->h/2)-roi->h/2;
  }
  return 0;
 
}
/*****************************************************************************/
 
/*****************************************************************************/
void roi_print(
  ROI *roi
) {
  int i;
 
  printf("imgfile: %s\n", roi->imgfile);
  printf("zoom: %f\n  recon_zoom: %f\n", roi->zoom, roi->recon_zoom);
#if 0
  printf("matnum: %d ; plane: %d frame: %d\n",
    roi->matnum, roi_mplane(roi->matnum), roi_mframe(roi->matnum) );
#endif
  printf("type: %d\n status: %d\n", roi->type, roi->status);
  printf("pos_x=%d pos_y=%d w=%d h=%d t=%d\n",
    roi->pos_x, roi->pos_y, roi->w, roi->h, roi->t);
  printf("roi=%d _%s_\n", roi->roi, roi->roiname);
  printf("points=%d\n", roi->point_nr);
  if(roi->point_nr>0) {
    for(i=0; i<roi->point_nr; i++) printf("(%d,%d) ", roi->x[i], roi->y[i]);
    printf("\n");
  }
}
/*****************************************************************************/
 
/*****************************************************************************/
typedef struct {
  float minY;
  float maxY;
  float x;
  float m;
} PolygonEdge;
 
static void fill_traceroi(
  char **m,
  int dimx,
  int dimy,
  int pos_x,
  int pos_y,
  int *roix,
  int *roiy,
  int point_nr
) {
  PolygonEdge *edges=malloc(sizeof(PolygonEdge)*point_nr);
  int edge=0;
  for(int p=0; p<point_nr; p++) {
    float maxX;
    int n;
    if(p==point_nr-1) n=0; else n=p+1;
    if(roiy[p]>roiy[n]) {
      edges[edge].minY=roiy[n]+pos_y;
      edges[edge].x=roix[n]+pos_x;
      edges[edge].maxY=roiy[p]+pos_y;
      maxX=roix[p]+pos_x;
    } else if(roiy[p]<roiy[n]) {
      edges[edge].minY=roiy[p]+pos_y;
      edges[edge].x=roix[p]+pos_x;
      edges[edge].maxY=roiy[n]+pos_y;
      maxX=roix[n]+pos_x;
    } else 
      continue;
    edges[edge].m=(maxX-edges[edge].x)/(edges[edge].maxY-edges[edge].minY);
    /* printf("edge[%d]: p=%d minY=%g maxY=%g x=%g m=%g\n", edge, p, edges[edge].minY, 
              edges[edge].maxY, edges[edge].x, edges[edge].m); */
    edge++;
  }
  int *px=malloc(sizeof(int)*point_nr);
  for(int y=0; y<dimy; y++) {
    //printf("y=%d\n", y);
    int pc=0;
    for(int p=0; p<edge; p++) {
//    if(temp_roundf(edges[p].minY)<=y && y<temp_roundf(edges[p].maxY)) {
//    2018-12-20 by VO
      if(temp_roundf(edges[p].minY)<=y && y<=temp_roundf(edges[p].maxY)) {
        px[pc]=temp_roundf(edges[p].x+edges[p].m*(y-temp_roundf(edges[p].minY)));
        if(px[pc]<0) px[pc]=0;
        pc++;
      }
    }
    int odd=0;
    for(int x=0; x<dimx; x++) {
      //printf(" x=%d odd=%d\n", x, odd);
      m[x][y]=odd;
      for(int p=0; p<pc; p++) {
        if(x==px[p]) {
          if(pc>1) odd=!odd;
          m[x][y]=1;
        }
      }
      //printf("  m[%d][%d]=%d\n", x, y, m[x][y]);
    }
  }
  free(px);
  free(edges);
}
 
int roi_onoff(
  ROI *roi,
  int dimx,
  int dimy,
  char **m
) {
  /* Check the parameters */
  if(roi==NULL || dimx<2 || dimy<2 || m==NULL) {
    strcpy(roierrmsg, "invalid arguments for roi_onoff()");
    return 1;
  }
  /* Prepare */
  float pos_x=roi->pos_x;
  float pos_y=roi->pos_y;
  for(int x=0; x<dimx; x++) for(int y=0; y<dimy; y++) m[x][y]=0;
 
  /* Fill the matrix */
  switch(roi->type) {
    case ROI_TRACE:    /* Trace ROI */
        fill_traceroi(m, dimx, dimy, pos_x, pos_y, roi->x, roi->y, roi->point_nr);
      break;
    case ROI_RECTANGULAR: /* Rectangle ROI */
      for(int x=pos_x; x<pos_x+roi->w; x++)
        for(int y=pos_y; y<pos_y+roi->h; y++)
          m[x][y]=1;
      break;
    case ROI_CIRCULAR:{    /* Circle ROI */
      int x1, x2;
      int y1, y2;
      int r=roi->w;
      for(int y=0; y<r; y++) {
        y1=pos_y-y;
        y2=pos_y+y;
        x1=temp_roundf(sqrt(r*r-y*y));
        x2=pos_x+x1;
        x1=pos_x-x1;
        if(x1<0) x1=0;
        if(x2>dimx) x2=dimx;
        while(x1<x2) {
          if(y1>=0 && y1<dimy) m[x1][y1]=1;
          if(y2>=0 && y2<dimy) m[x1][y2]=1;
          x1++;
        }
      }
    } break;
    case ROI_ELLIPSE: { /* Ellipse ROI */
      int x1, x2, a;
      int y1, y2, b;
      a=roi->w;
      b=roi->h;
      for(int y=0; y<b; y++) {
        y1=pos_y-y;
        y2=pos_y+y;
        x1=sqrt(1-(y*y)/(float)(b*b))*a;
        x2=pos_x+x1;
        x1=pos_x-x1;
        if(x1<0) x1=0;
        if(x2>dimx) x2=dimx;
        while(x1<x2) {
          if(y1>=0 && y1<dimy) m[x1][y1]=1;
          if(y2>=0 && y2<dimy) m[x1][y2]=1;
          x1++;
        }
      }
    } break;
    default:
      strcpy(roierrmsg,"invalid roi type");
      return 2;
  }
  return 0;
}
/*****************************************************************************/
int roi_save(
  const char *fname, 
  ROI_list *rl
) {
  RoiList *rois;
  FILE *fp;
 
  /* Check arguments */
  if(!fname[0] || rl==NULL || rl->nr<1) {
    strcpy(roierrmsg, "invalid arguments for roiSave()");
    return 1;
  }
 
  /* Create file */
  if((fp=fopen(fname, "w"))==NULL) {
    strcpy(roierrmsg, "cannot open file");
    return 2;
  }
  /* Loop and save the ROIs */
  rois=rl->rois;
  while(rois) {
    if(roi_append(fp,rois->roi)) return 3;
    rois=rois->next;
  }
  return 0;
}
/*****************************************************************************/
 
/*****************************************************************************/
int roi_append(
  FILE *fp,
  ROI *roi
) {
  /* Write first line */
  if(fprintf(fp, "*%s %f %f %d %d %d %d %d %d %d %d %d %s///%d %d\n",
    roi->imgfile, roi->zoom, roi->recon_zoom,
    roi->matnum, roi->type, roi->status,
    roi->pos_x, roi->pos_y, roi->w, roi->h,
    roi->t, roi->roi, roi->roiname,
    0, roi->point_nr ) < 10 ) 
  {
    strcpy(roierrmsg, "cannot write data");
    fclose(fp); return 1;
  }
  /* Write second line, if trace ROI */
  if(roi->type==ROI_TRACE) {
    int j;
    for(j=0; j<roi->point_nr; j++)
      fprintf(fp, "%d %d ", roi->x[j], roi->y[j]);
    fprintf(fp, "\n");
  }
  return 0;
}
 
/*****************************************************************************/
int roi_append_n(
  FILE *fp,
  ROI_list *rl,
  int ind
) {
  RoiList *rois=rl->rois;
  int pos=0;
  while(pos<ind) {
    rois=rois->next;
    pos++;
  }
  return roi_append(fp,rois->roi);
}
/*****************************************************************************/
 
/*****************************************************************************/
int jsqrt(
  int n
) {
  int x, x2, q, r, t;
 
  if(n<2) return n;
  t=x=n; while(t>>=2) x>>=1; x++;
  while(1) {
    q=n/x; r=n%x; 
    if(x<=q) {x2=x+2; if(q<x2 || (q==x2 && r==0)) break;}
    x=(x+q)>>1;
  }
  return x;
}
/*****************************************************************************/
 
/*****************************************************************************/
int roi_mplane(
  int matnum
) {
  Matval matval;
 
  mat_numdoc(matnum, &matval);
  return matval.plane;
}
/*****************************************************************************/
 
/*****************************************************************************/
int roi_mframe(
  int matnum
) {
  Matval matval;
 
  mat_numdoc(matnum, &matval);
  return matval.frame;
}
/*****************************************************************************/
 
/*****************************************************************************/