convo.c

卷积程序

A2F(InputStruct * In_Ptr, double **A_rz)
{
  short       nz = In_Ptr->nz, nr = In_Ptr->nr;
  double      dr = In_Ptr->dr, dz = In_Ptr->dz;
  short       ir, iz;
  double      mua;

  for (ir = 0; ir < nr; ir++)
    for (iz = 0; iz < nz; iz++) {
      mua = In_Ptr->layerspecs[IzToLayer(iz, In_Ptr)].mua;
      if (mua > 0.0)
	A_rz[ir][iz] /= mua;
    }
}

/****************************************************************
 *	Fluence to absorption density. F = A*mua;
 ****/
void 
F2A(InputStruct * In_Ptr, double **A_rz)
{
  short       nz = In_Ptr->nz, nr = In_Ptr->nr;
  double      dr = In_Ptr->dr, dz = In_Ptr->dz;
  short       ir, iz;
  double      mua;

  for (ir = 0; ir < nr; ir++)
    for (iz = 0; iz < nz; iz++) {
      mua = In_Ptr->layerspecs[IzToLayer(iz, In_Ptr)].mua;
      if (mua > 0.0)
	A_rz[ir][iz] *= mua;
    }
}

/****************************************************************
 ****/
void 
BranchContOrigCmd(char *Cmd_Str,
		  InputStruct * In_Ptr,
		  OutStruct * Out_Ptr)
{
  char        ch;

  switch (toupper(Cmd_Str[0])) {
  case 'A':
    IsoPlot(Out_Ptr->A_rz, In_Ptr->nr - 1, In_Ptr->nz - 1,
	    In_Ptr->dr, In_Ptr->dz);
    break;
  case 'F':
    A2F(In_Ptr, Out_Ptr->A_rz);
    IsoPlot(Out_Ptr->A_rz, In_Ptr->nr - 1, In_Ptr->nz - 1,
	    In_Ptr->dr, In_Ptr->dz);
    F2A(In_Ptr, Out_Ptr->A_rz);
    break;
  case 'R':
    IsoPlot(Out_Ptr->Rd_ra, In_Ptr->nr - 1, In_Ptr->na - 1,
	    In_Ptr->dr, In_Ptr->da);
    break;
  case 'T':
    IsoPlot(Out_Ptr->Tt_ra, In_Ptr->nr - 1, In_Ptr->na - 1,
	    In_Ptr->dr, In_Ptr->da);
    break;
  case 'H':
    ShowContOrigMenu(In_Ptr->in_fname);
    break;
  case 'Q':
    break;
  default:
    puts("...Wrong command");
  }
}
/****************************************************************
 ****/
void 
ContourOrigData(InputStruct * In_Ptr,
		OutStruct * Out_Ptr)
{
  char        cmd_str[STRLEN];

  if (!Out_Ptr->allocated)
    puts("...No data to output");
  else
    do {
      printf("\n> Contour output of mcml data (h for help) => ");
      do
	gets(cmd_str);
      while (!strlen(cmd_str));	/* avoid null string. */
      BranchContOrigCmd(cmd_str, In_Ptr, Out_Ptr);
    } while (toupper(cmd_str[0]) != 'Q');
}

/****************************Scanning***************************/
/****************************************************************
 ****/
void 
ShowScanOrigMenu(char *in_fname)
{
  printf("Ar = absorption vs r @ fixed z [1/cm3]\n");
  printf("Az = absorption vs z @ fixed r [1/cm3]\n");
  printf("Fr = fluence vs r @ fixed z [1/cm2]\n");
  printf("Fz = fluence vs z @ fixed r [1/cm2]\n");
  printf("Rr = diffuse reflectance vs r @ fixed angle [1/(cm2 sr)]\n");
  printf("Ra = diffuse reflectance vs angle @ fixed r [1/(cm2 sr)]\n");
  printf("Tr = transmittance vs r @ fixed angle [1/(cm2 sr)]\n");
  printf("Ta = transmittance vs angle @ fixed r [1/(cm2 sr)]\n");
  printf("Q  = quit\n");
  printf("* input filename: %s \n", in_fname);
}

/****************************************************************
 *	Ext is either "Ars" or "Frs".
 ****/
void 
ScanOrigA_r(char *Ext, InputStruct * In_Ptr, double **A_rz)
{
  short       ir, iz, nr = In_Ptr->nr, nz = In_Ptr->nz;
  double      r, z, dr = In_Ptr->dr, dz = In_Ptr->dz;
  FILE       *file;

  file = GetWriteFile(Ext);
  if (file == NULL)
    return;

  printf("z grid separation is %-10.4lg cm.\n", dz);
  printf("Input fixed z index (0 - %2hd): ", nz - 1);
  iz = GetShort(0, nz - 1);
  fprintf(file, "%-12s\t%-s@z=%-9.3lg\n", "r[cm]", Ext, dz * (iz + 0.5));
  for (ir = 0; ir < nr; ir++) {
    r = (ir + 0.5) * dr;
    fprintf(file, "%-12.4E\t%-12.4E\n", r, A_rz[ir][iz]);
  }

  fclose(file);
}

/****************************************************************
 *	Ext is either "Azs" or "Fzs".
 ****/
void 
ScanOrigA_z(char *Ext, InputStruct * In_Ptr, double **A_rz)
{
  short       ir, iz, nr = In_Ptr->nr, nz = In_Ptr->nz;
  double      r, z, dr = In_Ptr->dr, dz = In_Ptr->dz;
  FILE       *file;

  file = GetWriteFile(Ext);
  if (file == NULL)
    return;

  printf("r grid separation is %-10.4lg cm.\n", dr);
  printf("Input fixed r index (0 - %2hd): ", nr - 1);
  ir = GetShort(0, nr - 1);
  fprintf(file, "%-12s\t%-s@r=%-9.3lg\n", "z[cm]", Ext, dr * (ir + 0.5));
  for (iz = 0; iz < nz; iz++) {
    z = (iz + 0.5) * dz;
    fprintf(file, "%-12.4E\t%-12.4E\n", z, A_rz[ir][iz]);
  }

  fclose(file);
}

/****************************************************************
 ****/
void 
ScanOrigRd_r(InputStruct * In_Ptr, double **Rd_ra)
{
  short       ir, ia, nr = In_Ptr->nr, na = In_Ptr->na;
  double      r, a, dr = In_Ptr->dr, da = In_Ptr->da;
  FILE       *file;
  char        fname[STRLEN];

  strcpy(fname, "Rrs");
  file = GetWriteFile(fname);
  if (file == NULL)
    return;

  printf("Angle grid separation is %-10.4lg rad.\n", da);
  printf("Input fixed angle index (0 - %2hd): ", na - 1);
  ia = GetShort(0, na - 1);
  fprintf(file, "%-12s\t%-s@a=%-9.3lg\n", "r[cm]", fname, da * (ia + 0.5));
  for (ir = 0; ir < nr; ir++) {
    r = (ir + 0.5) * dr;
    fprintf(file, "%-12.4E\t%-12.4E\n", r, Rd_ra[ir][ia]);
  }

  fclose(file);
}

/****************************************************************
 ****/
void 
ScanOrigRd_a(InputStruct * In_Ptr, double **Rd_ra)
{
  short       ir, ia, nr = In_Ptr->nr, na = In_Ptr->na;
  double      r, a, dr = In_Ptr->dr, da = In_Ptr->da;
  FILE       *file;
  char        fname[STRLEN];

  strcpy(fname, "Ras");
  file = GetWriteFile(fname);
  if (file == NULL)
    return;

  printf("r grid separation is %-10.4lg cm.\n", dr);
  printf("Input fixed r index (0 - %2hd): ", nr - 1);
  ir = GetShort(0, nr - 1);
  fprintf(file, "%-12s\t%-s@r=%-9.3lg\n", "a[rad]", fname, dr * (ir + 0.5));
  for (ia = 0; ia < na; ia++) {
    a = (ia + 0.5) * da;
    fprintf(file, "%-12.4E\t%-12.4E\n", a, Rd_ra[ir][ia]);
  }

  fclose(file);
}

/****************************************************************
 ****/
void 
ScanOrigTt_r(InputStruct * In_Ptr, double **Tt_ra)
{
  short       ir, ia, nr = In_Ptr->nr, na = In_Ptr->na;
  double      r, a, dr = In_Ptr->dr, da = In_Ptr->da;
  FILE       *file;
  char        fname[STRLEN];

  strcpy(fname, "Trs");
  file = GetWriteFile(fname);
  if (file == NULL)
    return;

  printf("Angle grid separation is %-10.4lg rad.\n", da);
  printf("Input fixed angle index (0 - %2hd): ", na - 1);
  ia = GetShort(0, na - 1);
  fprintf(file, "%-12s\t%-s@a=%-9.3lg\n", "r[cm]", fname, da * (ia + 0.5));
  for (ir = 0; ir < nr; ir++) {
    r = (ir + 0.5) * dr;
    fprintf(file, "%-12.4E\t%-12.4E\n", r, Tt_ra[ir][ia]);
  }

  fclose(file);
}

/****************************************************************
 ****/
void 
ScanOrigTt_a(InputStruct * In_Ptr, double **Tt_ra)
{
  short       ir, ia, nr = In_Ptr->nr, na = In_Ptr->na;
  double      r, a, dr = In_Ptr->dr, da = In_Ptr->da;
  FILE       *file;
  char        fname[STRLEN];

  strcpy(fname, "Tas");
  file = GetWriteFile(fname);
  if (file == NULL)
    return;

  printf("r grid separation is %-10.4lg cm.\n", dr);
  printf("Input fixed r index (0 - %2hd): ", nr - 1);
  ir = GetShort(0, nr - 1);
  fprintf(file, "%-12s\t%-s@r=%-9.3lg\n", "a[rad]", fname, dr * (ir + 0.5));
  for (ia = 0; ia < na; ia++) {
    a = (ia + 0.5) * da;
    fprintf(file, "%-12.4E\t%-12.4E\n", a, Tt_ra[ir][ia]);
  }

  fclose(file);
}

/****************************************************************
 ****/
void 
BranchScanOrigA(char *Cmd_Str,
		InputStruct * In_Ptr,
		OutStruct * Out_Ptr)
{
  char        fname[STRLEN];

  switch (toupper(Cmd_Str[1])) {
  case 'R':
    strcpy(fname, "Ars");
    ScanOrigA_r(fname, In_Ptr, Out_Ptr->A_rz);
    break;
  case 'Z':
    strcpy(fname, "Azs");
    ScanOrigA_z(fname, In_Ptr, Out_Ptr->A_rz);
    break;
  default:
    puts("...Wrong command");
  }
}

/****************************************************************
 ****/
void 
BranchScanOrigF(char *Cmd_Str,
		InputStruct * In_Ptr,
		OutStruct * Out_Ptr)
{
  char        fname[STRLEN];

  A2F(In_Ptr, Out_Ptr->A_rz);

  switch (toupper(Cmd_Str[1])) {
  case 'R':
    strcpy(fname, "Frs");
    ScanOrigA_r(fname, In_Ptr, Out_Ptr->A_rz);
    break;
  case 'Z':
    strcpy(fname, "Fzs");
    ScanOrigA_z(fname, In_Ptr, Out_Ptr->A_rz);
    break;
  default:
    puts("...Wrong command");
  }

  F2A(In_Ptr, Out_Ptr->A_rz);
}

/****************************************************************
 ****/
void 
BranchScanOrigR(char *Cmd_Str,
		InputStruct * In_Ptr,
		OutStruct * Out_Ptr)
{
  switch (toupper(Cmd_Str[1])) {
    case 'R':
    ScanOrigRd_r(In_Ptr, Out_Ptr->Rd_ra);
    break;
  case 'A':
    ScanOrigRd_a(In_Ptr, Out_Ptr->Rd_ra);
    break;
  default:
    puts("...Wrong command");
  }
}

/****************************************************************
 ****/
void 
BranchScanOrigT(char *Cmd_Str,
		InputStruct * In_Ptr,
		OutStruct * Out_Ptr)
{
  switch (toupper(Cmd_Str[1])) {
    case 'R':
    ScanOrigTt_r(In_Ptr, Out_Ptr->Tt_ra);
    break;
  case 'A':
    ScanOrigTt_a(In_Ptr, Out_Ptr->Tt_ra);
    break;
  default:
    puts("...Wrong command");
  }
}

/****************************************************************
 ****/
void 
BranchScanOrigCmd(char *Cmd_Str,
		  InputStruct * In_Ptr,
		  OutStruct * Out_Ptr)
{
  char        ch;

  switch (toupper(Cmd_Str[0])) {
  case 'A':
    BranchScanOrigA(Cmd_Str, In_Ptr, Out_Ptr);
    break;
  case 'F':
    BranchScanOrigF(Cmd_Str, In_Ptr, Out_Ptr);
    break;
  case 'R':
    BranchScanOrigR(Cmd_Str, In_Ptr, Out_Ptr);
    break;
  case 'T':
    BranchScanOrigT(Cmd_Str, In_Ptr, Out_Ptr);
    break;
  case 'H':
    ShowScanOrigMenu(In_Ptr->in_fname);
    break;
  case 'Q':
    break;
  default:
    puts("...Wrong command");
  }
}

/****************************************************************
 ****/
void 
ScanOrigData(InputStruct * In_Ptr,
	     OutStruct * Out_Ptr)
{
  char        cmd_str[STRLEN];

  if (!Out_Ptr->allocated)
    puts("...No data to output");
  else
    do {
      printf("\n> Scans of mcml data (h for help) => ");
      do
	gets(cmd_str);
      while (!strlen(cmd_str));	/* avoid null string. */
      BranchScanOrigCmd(cmd_str, In_Ptr, Out_Ptr);
    } while (toupper(cmd_str[0]) != 'Q');
}