The analysis macro, scan6.mac, MUST BE EDITED SLIGHTLY PRIOR TO USE!!! During failed attempts to load this macro onto the web I discovered that deletions to the code would occur in all statements where any letter followed a "less than" sign (<). This occurs because of incompatibility of coding statements using this symbol with its use for formatting statements in html. Therefore, in the code for this macro, any line that contained a "less than" sign followed by a letter, the sign was replaced with the words "less than". Users of these macros need to download the text of the code, print out a hard copy and find all such statements using the words "less than" and replace them with the "less than" symbol.

A general note about both macros: Since they were generated using OPTIMAS 3.10, and probably no one uses this version of the program any longer, it is likely that users with more recent versions of OPTIMAS will have to debug these macros and edit them to make them compatible with the updated versions of ALI.


Scan6.mac - Analysis macro


DEFINE Find_polyp_box (vPolyp_box, dx_step, rDist)
{
  // The ROI should be right when coming to this point!!
  // Isolate and separate polyp base from edge of ROI
  // by putting a black line at the right edge, and along margin of ROI.
  PutPixelRect(vPix_r, vBlackline);
  PutPixelRect(vPix_l, vBlackline);
  PutPixelRect(vPix_t, vBlackline);
  PutPixelRect(vPix_b, vBlackline);

  //Set area sampling parameters, delete areas touching ROI and nested areas
  AreaCNVFactors=0.0 : 1.0 : -1.0 : 256.0 : 0.0 : 100.0 : 1.0 : 1.0 : -1.0;
  CreateArea( , FALSE, TRUE);
  Classify( );
  ObjectClass( Class1, TRUE); //edit selected areas so that only right-sized objects included
  SetExport(ArSampledPoints, 0, TRUE);

  MultipleMode = TRUE;
  DataCollectionType = 2;
  MultipleExtract( );
  SetExport(ArSampledPoints, 0, FALSE); 

// Place selected area objects in frame buffer
ROIToList( , "ARM_Original");

// Copy area screen object interiors into a bitmap (fills in polyps so they are all white)
ARM_BitMap = BitMapCreate( , 0, 255, TRUE, 0);

// Return this newly created bitmap to the display which now has an all black background
BitmapProcess(ARM_BitMap, , 22);
BitmapDelete(ARM_BitMap);

DeleteImage("ARM_Original");
ObjectWildCardList("ARM_*", 2);

// Dilate the image to enclose broken margins into a complete polyp outline
InvertFilter( );
DilateFilter( );


  // Isolate and separate polyp base from edge of ROI
  // by putting a black line at the right edge, and along margin of ROI.
  PutPixelRect(vPix_r, vBlackline);
  PutPixelRect(vPix_l, vBlackline);
  PutPixelRect(vPix_t, vBlackline);
  PutPixelRect(vPix_b, vBlackline);

  CreateArea( , FALSE, TRUE);
  Classify( );
  SetExport(ArSampledPoints, 0, TRUE);
  SetExport(ArArea, 1, TRUE);
  SetExport(ArPerimeter, 1, TRUE);

  Extract( );

  SetExport(ArSampledPoints, 0, FALSE);
  SetExport(ArArea, 1, FALSE);
  SetExport(ArPerimeter, 1, FALSE);

  // We use here as the unit 1000*(um)^2 which is _not_ (mm)^2, be careful!!!!
  rDist[1]=ArArea*1.849;
  // We use here as the unit 10um which is 10^(-5)m.  The other units are um!!!
  rDist[2]=ArPerimeter*4.3;

REAL vPoints = ArSampledPoints;
  vPolyp_box[0,0]=Min(vPoints[ ,0]);
  vPolyp_box[0,1]=Max(vPoints[ ,1]) + dy;
  vPolyp_box[1,0]=Max(vPoints[ ,0]);
  vPolyp_box[1,1]=Min(vPoints[ ,1]) - dy;

  dx_step=(vPolyp_box[1,0]-vPolyp_box[0,0])/iN_dx_step;

}; //end of function Find polyp box

DEFINE Ver_polyp_scan (vPolyp_box, iN_width_meas, iN_max_width_meas, dx_step, vPoints_polyp_axis,vDir_n)
{ //beginning of function
  REAL vLine[2,2]=(Max(vPolyp_box[ ,0])-dx_step*0.5) : Min(vPolyp_box[ ,1]) :: (Max(vPolyp_box[ ,0])-dx_step*0.5) : Max(vPolyp_box[ ,1]);
  do {
    // establish direction for successive line scans to proceed
    if (vDir_n[0]<0) vLine_ver_scan=vLine - (dx_step : 0.0)*iN_width_meas;
    else vLine_ver_scan=vLine + (dx_step : 0.0)*iN_width_meas;

    // Create a line and sample its luminance
    CreateLine(vLine_ver_scan);
    SetExport(LnLuminance, 0, TRUE);
    Extract( );
    SetExport(LnLuminance, 0, FALSE);
    // scan through this line pixel by pixel, from bottom to top, until
    // total line luminance is reached, pixel is extracted repeat process to get pixel on bottom margin,
    // midpoint y-coordinate along the polyp axis is obtained by geometric midpoint between luminance margins

     INTEGER iNs = 2*VectorLength(LnLuminance);
      lSumlum = Sum(LnLuminance);

      //scan through the line to determine edges of the polyp by
      //determining the coordinates at which all the luminance of the
      //vLine_P_Diam has been accounted for 
      INTEGER iN_pixel=0;
      lSum=0;
      do
	lSum = lSum + LnLuminance[iN_pixel++];
      while (lSum less than lSumlum);
      iN_pixel--;

      //scan through the line in the other direction performing the same
      //operation
      INTEGER iN_pix=iN_pixel;
      lSumlum = Sum(LnLuminance);
      lSum = 0;
      do
	lSum = lSum + LnLuminance[iN_pix--];
      while (lSum less than lSumlum);
      iN_pix++;
	 REAL y_on_axis=vLine_ver_scan[0,1] + ((REAL)iN_pixel+(REAL)iN_pix)/((REAL)iNs)*(vLine_ver_scan[1,1]-vLine_ver_scan[0,1]);
    // Both the x and y midpoint coordinates along the polyp axis encoded
    // in the vPoints_polyp_axis vector
    vPoints_polyp_axis[iN_width_meas, ]=vLine_ver_scan[0,0] : y_on_axis;
    ClearScreen(FALSE);
    //For second, and subsequent vertical scans along polyp, (1) create
    //line containing midpoint axis coordinates, and (2) establish
    //perpendicular line to polyp axis to measure diameter at each interval
    //along the axis where a vertical scan is done
    if (iN_width_meas>=1) {
      REAL vDir=vPoints_polyp_axis[iN_width_meas, ] - vPoints_polyp_axis[iN_width_meas-1, ];
      vDir_n=vDir/Sqrt(vDir[0]*vDir[0]+vDir[1]*vDir[1]);
      //normalized directional vector of polyp axis at that width measure

      REAL vPerp[2] = -vDir[1] : vDir[0];
      //perpendicular directional vector to the polyp axis at that width
      //measure scale factors that determine how far (given the slope of
      //the perpendicular) the top and bottom of the polyp box is from the
      //current midpoints of the polyp axis
      REAL x_low = vPolyp_box[0,0];
      REAL x_high = vPolyp_box[1,0];
      REAL y_low = vPolyp_box[1,1];
      REAL y_high = vPolyp_box[0,1];
      REAL mid_x = 0.5*(vPoints_polyp_axis[iN_width_meas,0]+vPoints_polyp_axis[iN_width_meas-1,0]);
      REAL mid_y = 0.5*(vPoints_polyp_axis[iN_width_meas,1]+vPoints_polyp_axis[iN_width_meas-1,1]);
      REAL r=0;
      // First try hitting the higher line of the box
      // y_high = r*vPerp[1] + mid_y
      REAL hit_high_y=0;
      r=(y_high-mid_y)/vPerp[1];
      // if that is outside of the box, recalculate
      hit_high_x = mid_x + r * vPerp[0];
      if (hit_high_x > x_high) r=(x_high-mid_x)/vPerp[0];
      else if (hit_high_x < x_low) r=(x_low-mid_x)/vPerp[0];
      hit_high_y = mid_y + r * vPerp[1];
      hit_high_x = mid_x + r * vPerp[0];
      // Then try hitting the lower line of the box
      // y_low = r*vPerp[1] + mid_y
      REAL hit_low_y=0;
      r=(y_low-mid_y)/vPerp[1];
      // if that is outside of the box, recalculate
      hit_low_x = mid_x + r * vPerp[0];
      if (hit_low_x > x_high) r=(x_high-mid_x)/vPerp[0];
      else if (hit_low_x < x_low) r=(x_low-mid_x)/vPerp[0];
      hit_low_y = mid_y + r * vPerp[1];
      hit_low_x = mid_x + r * vPerp[0];

      //Create line given the midpoint of the axis at that width measure,
      //the perpendicular slope and the top and bottom scale factors that
      //give the x & y coords of the endpoints of the line along the
      //margins of the polyp box 
      REAL vLine_P_Diam[2,2] = hit_high_x : hit_high_y :: hit_low_x : hit_low_y;
      CreateLine(vLine_P_Diam);
      SetExport(LnLuminance, 0, TRUE);
      SetExport(LnLength, 0, TRUE);
      Extract( );
      SetExport(LnLuminance, 0, FALSE);
      SetExport(LnLength, 0, FALSE);

      //Extract the luminance profile from line and determine the vector
      //length in pixels (how many luminance values calculated that fill
      //the vector corresponding to this line)
      INTEGER iNs = VectorLength(LnLuminance);
      lSumlum = Sum(LnLuminance);

      //scan through the line to determine edges of the polyp by
      //determining the coordinates at which all the luminance of the
      //vLine_P_Diam has been accounted for 
      INTEGER iN_pixel=0;
      lSum=0;
      do
	lSum = lSum + LnLuminance[iN_pixel++];
      while (lSum less than lSumlum);
      iN_pixel--;

      //scan through the line in the other direction performing the same
      //operation
      INTEGER iN_pix=iN_pixel;
      lSumlum = Sum(LnLuminance);
      lSum = 0;
      do
	lSum = lSum + LnLuminance[iN_pix--];
      while (lSum less than lSumlum);
      iN_pix++;
      REAL Frac = 1.0*(iN_pixel - iN_pix)/iNs;
      // We have 1204um and _not_ 28cm horizontally!!!
      vLine_polyp_diam[iN_width_meas-1] = LnLength*Frac*43.0;
    }
  } while(Abs(vDir_n[0])>0.6 && ++iN_width_meas less than iN_max_width_meas);
}  //end of vertical scan function

hLib=LoadMacroLibrary("matrix.oml");

DEFINE Smooth_Points_polyp_axis (iOrder, vSmooth_Points_polyp_axis)
{ //beginning of function
  hLib=LoadMacroLibrary("matrix.oml");
  REAL vXval=vPoints_polyp_axis[, 0];
  REAL vYval=vPoints_polyp_axis[, 1];
  REAL vCoef=FitToPolynomial(vXval, vYval, iOrder);
  REAL vYval_smooth=ExpandPolynomial(vXval, vCoef);
  vSmooth_Points_polyp_axis[, 0] = vXval;
  vSmooth_Points_polyp_axis[, 1] = vYval_smooth;
  CreateLine(vSmooth_Points_polyp_axis);
} //end of function

DEFINE Find_Polyp_tip (vPoints_polyp_axis, vSmooth_points_polyp_axis, dx, dy, vPolyp_box, iN_width_meas, lSumlum, vDir_n)
{  //beginning of function
  REAL vDir=vSmooth_points_polyp_axis[iN_width_meas-1, ] - vSmooth_Points_polyp_axis[iN_width_meas-2,];
  vDir_n=vDir/Sqrt(vDir*vDir);
  REAL vLine[2,2];
  if (ABS(vDir_n[0])<=0.7) {
    vLine = Min(vPolyp_box[ ,0]):vSmooth_Points_polyp_axis[iN_width_meas-1,1] :: Max(vPolyp_box[ ,0]):vSmooth_Points_polyp_axis[iN_width_meas-1, 1];

    // begin horizontal scan, first check direction to proceed
    INTEGER i=0;
    REAL vLine_hor_scan[2,2];
    do {
      if (vDir_n[1]<0) vLine_hor_scan=vLine - (0.0 : dy)*i++;
      else vLine_hor_scan=vLine + (0.0 : dy)*i++;
      CreateLine(vLine_hor_scan);
      SetExport(LnLuminance, 0, TRUE);
      Extract( );
      SetExport(LnLuminance, 0, FALSE);
      lSumlum=Sum(LnLuminance);
    } while (lSumlum>250);  //scan line out of polyp body

    if (vDir_n[1]<0) vLine_hor_scan=vLine - (0.0:dy)*(i-2);
    else vLine_hor_scan=vLine + (0.0:dy)*(i-2);
    CreateLine(vLine_hor_scan);
    SetExport(LnLuminance, 0, TRUE);
    Extract( );
    SetExport(LnLuminance, 0, FALSE);
    lSumlum=Sum(LnLuminance)/2;
    LONG lSum=0;
    INTEGER iN_pixel=0;
    do
      lSum = lSum + LnLuminance[iN_pixel++];
    while (lSum less than lSumlum);
    REAL x_on_axis=vLine_hor_scan[0,0] + ((REAL)iN_pixel - 0.5)/((REAL)VectorLength(LnLuminance))*(vLine_ver_scan[1,0]-vLine_ver_scan[0,0]);
    vSmooth_Points_polyp_axis[iN_width_meas, ] = x_on_axis : vLine_hor_scan[0,1];
    // end of horizontal scan
  } else {
    // Do a vertical scan
    vLine=vSmooth_Points_polyp_axis[iN_width_meas-1, 0]:Min(vPolyp_box[ ,1]) ::vSmooth_Points_polyp_axis[iN_width_meas-1,0]:Max(vPolyp_box[ ,1]);
    //continue with vertical scan, first check direction to proceed
    INTEGER i=0;
    REAL vLine_ver_scan[2,2];
    do {
      if (vDir_n[0]<0) vLine_ver_scan = vLine - (dx : 0.0)*i++;
      else vLine_ver_scan = vLine + (dx : 0.0)*i++;
      CreateLine(vLine_ver_scan);
      SetExport(LnLuminance, 0, TRUE);
      Extract( );
      SetExport(LnLuminance, 0, FALSE);
      lSumlum = Sum(LnLuminance);
    } while (lSumlum > 250);  // scan line out of polyp body

    if (vDir_n[0]<0) vLine_ver_scan = vLine - (dx : 0.0)*(i-2);
    else vLine_ver_scan = vLine - (dx : 0.0)*(i-2);
    CreateLine(vLine_ver_scan);
    SetExport(LnLuminance, 0, TRUE);
    Extract( );
    SetExport(LnLuminance, 0, FALSE);
    lSumlum = Sum(LnLuminance)/2;
    INTEGER iN_pixel = 0;
    LONG lSum = 0;
    do 
      lSum = lSum + LnLuminance[iN_pixel++];
    while (lSum less than lSumlum);
    REAL y_on_axis = vLine_ver_scan[0,1] + ((REAL)iN_pixel - 0.5)/((REAL)VectorLength(LnLuminance))*(vLine_ver_scan[1,1]-vLine_ver_scan[0,1]);
    vSmooth_Points_polyp_axis[iN_width_meas, ] = vLine_ver_scan[0,0] : y_on_axis;
  } //end vertical scan for tip
} //end of function to find polyp tip
   
DEFINE Write_measurements_to_file (meas_file, rDist)
{ //beg of function
  df = OpenFile(meas_file, 0x0002);
  PositionFile(df, 0L, 2);
  CHAR Format = "%8.4f"::""::" "::"\r\n"::"";
  cDist = ToText(rDist, Format);
  WriteFile(df, cDist);
  CloseFile(df);
}  //end of function

DEFINE Polyp_axis_tracing ( )
{ // begin measurement block when image is loaded and binarized
  Find_Polyp_box (vPolyp_box, dx_step, rDist);
  // - - - - - - - - - - - - - - Set for Vertical Polyp scan
  INTEGER iN_width_meas = 0;
  GLOBAL LONG lSumlum = 0;
  GLOBAL REAL vDir_n[2] = -1 : 0;
  Ver_polyp_scan(vPolyp_box, iN_width_meas, iN_max_width_meas, dx_step,vPoints_polyp_axis, vDir_n );
  // - - - - - - - - - - - - - - Smooth points along polyp midpoint axis
  iOrder = 4;
  Smooth_Points_polyp_axis(iOrder, vSmooth_Points_polyp_axis );

  Find_Polyp_tip(vPoints_polyp_axis, vSmooth_points_polyp_axis, dx, dy, vPolyp_box, iN_max_width_meas, lSumlum, vDir_n);
  CreateLine(vSmooth_Points_polyp_axis);
  SetExport(LnLength, 1, TRUE);
  Extract( );
  SetExport(LnLength, 1, FALSE);
  rDist[0] = (REAL) iN;
  // The machine uses 28cm frames (horizontally, but we have 1204um, therefore:
  rDist[3] = 43.0*LnLength;
  for (i=0; i less than iN_max_width_meas; i++) rDist[i+4]=vLine_polyp_diam[i];

  Write_measurements_to_file(meas_file, rDist);
} //end of function

//end of all functions

Control MACRO
Gastrovascular dynamics
Dudgeon Lab Home
E-mail
CSUN Biology Home


Copyright 1999, Steve Dudgeon
Page last modified March 25, 1999