C++ Cube::close方法代码示例

// Main program
void IsisMain() {
  ProcessImportPds p;
  Pvl pdsLabel;
  UserInterface &ui = Application::GetUserInterface();
  FileName inFile = ui.GetFileName("FROM");

  p.SetPdsFile(inFile.expanded(), "", pdsLabel);

  QString filename = FileName(ui.GetFileName("FROM")).baseName();
  FileName toFile = ui.GetFileName("TO");
  
  apollo = new Apollo(filename);

  utcTime = (QString)pdsLabel["START_TIME"];  
  
  // Setup the output cube attributes for a 16-bit unsigned tiff
  Isis::CubeAttributeOutput cao;
  cao.setPixelType(Isis::Real);
  p.SetOutputCube(toFile.expanded(), cao);

  // Import image
  p.StartProcess();
  p.EndProcess();

  cube.open(toFile.expanded(), "rw");
  
  // Once the image is imported, we need to find and decrypt the code
  if (apollo->IsMetric() && FindCode())
    TranslateCode();
    
  CalculateTransform();
  // Once we have decrypted the code, we need to populate the image labels
  TranslateApolloLabels(filename, &cube);
  cube.close();
}

void IsisMain() {
  UserInterface &ui = Application::GetUserInterface();
  try {
    // Open the cube
    Cube cube;
    cube.open(ui.GetFileName("FROM"), "rw");

    //check for existing polygon, if exists delete it
    if(cube.label()->hasObject("Polygon")) {
      cube.label()->deleteObject("Polygon");
    }

    // Get the camera, interpolate to a parabola
    Camera *cam = cube.camera();
    if(cam->DetectorMap()->LineRate() == 0.0) {
      QString msg = "[" + ui.GetFileName("FROM") + "] is not a line scan camera";
      throw IException(IException::User, msg, _FILEINFO_);
    }
    cam->instrumentRotation()->SetPolynomial();

    // Get the instrument pointing keyword from the kernels group and update
    // its value to table.
    Isis::PvlGroup kernels =
      cube.label()->findGroup("Kernels", Isis::Pvl::Traverse);

    // Save original kernels in keyword before changing to "Table" in the kernels group
    PvlKeyword origCk = kernels["InstrumentPointing"];

    // Write out the "Table" label to the tabled kernels in the kernels group
    kernels["InstrumentPointing"] = "Table";

    // And finally write out the original kernels after Table
    for (int i = 0;  i < origCk.size();  i++) {
      kernels["InstrumentPointing"].addValue(origCk[i]);
    }

    cube.putGroup(kernels);

    // Pull out the pointing cache as a table and write it
    Table cmatrix = cam->instrumentRotation()->Cache("InstrumentPointing");
    cmatrix.Label().addComment("Smoothed using spicefit");
    cube.write(cmatrix);
    cube.close();
  }
  catch(IException &e) {
    QString msg = "Unable to fit pointing for [" + ui.GetFileName("FROM") + "]";
    throw IException(IException::User, msg, _FILEINFO_);
  }
}

void IsisMain() {
  // Create a serial number list
  UserInterface &ui = Application::GetUserInterface();
  QString filename = ui.GetFileName("FROM");
  SerialNumberList serialNumberList;
  serialNumberList.Add(filename);

  // Get the coordinate for updating the camera pointing
  // We will want to make the camera pointing match the lat/lon at this
  // line sample
  double samp1 = ui.GetDouble("SAMP1");
  double line1 = ui.GetDouble("LINE1");
  Latitude lat1(ui.GetDouble("LAT1"), Angle::Degrees);
  Longitude lon1(ui.GetDouble("LON1"), Angle::Degrees);
  Distance rad1;
  if(ui.WasEntered("RAD1")) {
    rad1 = Distance(ui.GetDouble("RAD1"), Distance::Meters);
  }
  else {
    rad1 = GetRadius(ui.GetFileName("FROM"), lat1, lon1);
  }

  // In order to use the bundle adjustment class we will need a control
  // network
  ControlMeasure * m = new ControlMeasure;
  m->SetCubeSerialNumber(serialNumberList.SerialNumber(0));
  m->SetCoordinate(samp1, line1);
//   m->SetType(ControlMeasure::Manual);
  m->SetType(ControlMeasure::RegisteredPixel);

  ControlPoint * p = new ControlPoint;
  p->SetAprioriSurfacePoint(SurfacePoint(lat1, lon1, rad1));
  p->SetId("Point1");
  p->SetType(ControlPoint::Fixed);
  p->Add(m);

  ControlNet cnet;
//  cnet.SetType(ControlNet::ImageToGround);
  cnet.AddPoint(p);

    // We need the target body
    Cube c;
    c.open(filename, "rw");
    //check for target name
    if(c.label()->hasKeyword("TargetName", PvlObject::Traverse)) {
//       c.Label()->findKeyword("TargetName");
      PvlGroup inst = c.label()->findGroup("Instrument", PvlObject::Traverse);
      QString targetName = inst["TargetName"];
      cnet.SetTarget(targetName);
    }
    c.close();

  // See if they wanted to solve for twist
  if(ui.GetBoolean("TWIST")) {
    double samp2 = ui.GetDouble("SAMP2");
    double line2 = ui.GetDouble("LINE2");
    Latitude lat2(ui.GetDouble("LAT2"), Angle::Degrees);
    Longitude lon2(ui.GetDouble("LON2"), Angle::Degrees);
    Distance rad2;
    if(ui.WasEntered("RAD2")) {
      rad2 = Distance(ui.GetDouble("RAD2"), Distance::Meters);
    }
    else {
      rad2 = GetRadius(ui.GetFileName("FROM"), lat2, lon2);
    }

    ControlMeasure * m = new ControlMeasure;
    m->SetCubeSerialNumber(serialNumberList.SerialNumber(0));
    m->SetCoordinate(samp2, line2);
    m->SetType(ControlMeasure::Manual);

    ControlPoint * p = new ControlPoint;
    p->SetAprioriSurfacePoint(SurfacePoint(lat2, lon2, rad2));
    p->SetId("Point2");
    p->SetType(ControlPoint::Fixed);
    p->Add(m);

    cnet.AddPoint(p);
  }

  // Bundle adjust to solve for new pointing
  try {
    BundleAdjust b(cnet, serialNumberList);
    b.SetSolveTwist(ui.GetBoolean("TWIST"));
    //    double tol = ui.GetDouble("TOL");
    //int maxIterations = ui.GetInteger("MAXITS");
    //b.Solve(tol, maxIterations);
    b.SetSolveCmatrix(BundleAdjust::AnglesOnly);
    b.SetSolveSpacecraftPosition(BundleAdjust::Nothing);
    b.SetErrorPropagation(false);
    b.SetOutlierRejection(false);
    b.SetSolutionMethod("SPECIALK");
    b.SetStandardOutput(true);
    b.SetCSVOutput(false);
    b.SetResidualOutput(true);
    b.SetConvergenceThreshold(ui.GetDouble("SIGMA0"));
    b.SetMaxIterations(ui.GetInteger("MAXITS"));

    b.SetDecompositionMethod(BundleAdjust::SPECIALK);
    b.SolveCholesky();
//.........这里部分代码省略.........

//.........这里部分代码省略.........
        QString testStartTime = (QString)timegrp["StartTime"];
        QString testStopTime = (QString)timegrp["StopTime"];
        if(testStartTime < startTime) {
          startTime = testStartTime;
          startClock = (QString)timegrp["SpacecraftClockStartCount"];
        }
        if(testStopTime > stopTime) {
          stopTime = testStopTime;
          stopClock = (QString)timegrp["spacecraftClockStopCount"];
        }
      }
    }

    //  Concatenate all TDI's and summing and specialProcessingFlat into one keyword
    PvlKeyword cpmmTdiFlag("cpmmTdiFlag");
    PvlKeyword cpmmSummingFlag("cpmmSummingFlag");
    PvlKeyword specialProcessingFlag("SpecialProcessingFlag");
    for(int i = 0; i < 14; i++) {
      cpmmTdiFlag += (QString)"";
      cpmmSummingFlag += (QString)"";
      specialProcessingFlag += (QString)"";
    }

    for(int i = 0; i < (int)clist.size(); i++) {
      Pvl *clab = clist[i]->label();
      PvlGroup cInst = clab->findGroup("Instrument", Pvl::Traverse);
      OriginalLabel cOrgLab;
      clist[i]->read(cOrgLab);
      PvlGroup cGrp = cOrgLab.ReturnLabels().findGroup("INSTRUMENT_SETTING_PARAMETERS", Pvl::Traverse);
      cpmmTdiFlag[(int)cInst["CpmmNumber"]] = (QString) cGrp["MRO:TDI"];
      cpmmSummingFlag[(int)cInst["CpmmNumber"]] = (QString) cGrp["MRO:BINNING"];

      if(cInst.hasKeyword("Special_Processing_Flag")) {
        specialProcessingFlag[cInst["CpmmNumber"]] = (QString) cInst["Special_Processing_Flag"];
      }
      else {
        // there may not be the keyword Special_Processing_Flag if no
        //keyword then set the output to NOMINAL
        specialProcessingFlag[cInst["CpmmNumber"]] = "NOMINAL";
      }
    }


    // Get the blob of original labels from first image in list
    OriginalLabel org;
    clist[0]->read(org);

    //close all cubes
    for(int i = 0; i < (int)clist.size(); i++) {
      clist[i]->close();
      delete clist[i];
    }
    clist.clear();

    // automos step
    QString list = ui.GetFileName("FROMLIST");
    QString toMosaic = ui.GetFileName("TO");
    QString MosaicPriority = ui.GetString("PRIORITY");

    QString parameters = "FROMLIST=" + list + " MOSAIC=" + toMosaic + " PRIORITY=" + MosaicPriority;
    ProgramLauncher::RunIsisProgram("automos", parameters);

    // write out new information to new group mosaic

    PvlGroup mos("Mosaic");
    mos += PvlKeyword("ProductId ", ProdId);
    mos += PvlKeyword(sourceProductId);
    mos += PvlKeyword("StartTime ", startTime);
    mos += PvlKeyword("SpacecraftClockStartCount ", startClock);
    mos += PvlKeyword("StopTime ", stopTime);
    mos += PvlKeyword("SpacecraftClockStopCount ", stopClock);
    mos += PvlKeyword("IncidenceAngle ", toString(Cincid), "DEG");
    mos += PvlKeyword("EmissionAngle ", toString(Cemiss), "DEG");
    mos += PvlKeyword("PhaseAngle ", toString(Cphase), "DEG");
    mos += PvlKeyword("LocalTime ", toString(ClocalSolTime), "LOCALDAY/24");
    mos += PvlKeyword("SolarLongitude ", toString(CsolarLong), "DEG");
    mos += PvlKeyword("SubSolarAzimuth ", toString(CsunAzimuth), "DEG");
    mos += PvlKeyword("NorthAzimuth ", toString(CnorthAzimuth), "DEG");
    mos += cpmmTdiFlag;
    mos += cpmmSummingFlag;
    mos += specialProcessingFlag;

    Cube mosCube;
    mosCube.open(ui.GetFileName("TO"), "rw");
    PvlObject &lab = mosCube.label()->findObject("IsisCube");
    lab.addGroup(mos);
    //add orginal label blob to the output cube
    mosCube.write(org);
    mosCube.close();

  }
  catch(IException &e) {
    for(int i = 0; i < (int)clist.size(); i++) {
      clist[i]->close();
      delete clist[i];
    }
    QString msg = "The mosaic [" + ui.GetFileName("TO") + "] was NOT created";
    throw IException(IException::User, msg, _FILEINFO_);
  }
} // end of isis main

void IsisMain() {
  // Get the list of cubes to stack
  Process p;
  UserInterface &ui = Application::GetUserInterface();
  FileList cubeList(ui.GetFileName("FROMLIST"));

  // Loop through the list
  int nsamps(0), nlines(0), nbands(0);
  PvlGroup outBandBin("BandBin");
  try {
    for(int i = 0; i < cubeList.size(); i++) {
      Cube cube;
      CubeAttributeInput inatt(cubeList[i].original());
      vector<QString> bands = inatt.bands();
      cube.setVirtualBands(bands);
      cube.open(cubeList[i].toString());
      if(i == 0) {
        nsamps = cube.sampleCount();
        nlines = cube.lineCount();
        nbands = cube.bandCount();
      }
      else {
        // Make sure they are all the same size
        if((nsamps != cube.sampleCount()) || (nlines != cube.lineCount())) {
          QString msg = "Spatial dimensions of cube [" +
                        cubeList[i].toString() + "] does not match other cubes in list";
          throw IException(IException::User, msg, _FILEINFO_);
        }
        // Get the total number of bands
        nbands += cube.bandCount();
      }

      // Build up the band bin group
      PvlObject &isiscube = cube.label()->findObject("IsisCube");
      if(isiscube.hasGroup("BandBin")) {
        PvlGroup &inBandBin = isiscube.findGroup("BandBin");
        for(int key = 0; key < inBandBin.keywords(); key++) {
          PvlKeyword &inKey = inBandBin[key];
          if(!outBandBin.hasKeyword(inKey.name())) {
            outBandBin += inKey;
          }
          else {
            PvlKeyword &outKey = outBandBin[inKey.name()];
            for(int index = 0; index < (int)inKey.size(); index++) {
              outKey.addValue(inKey[index], inKey.unit(index));
            }
          }
        }
      }
      cube.close();
    }
  }
  catch(IException &e) {
    QString msg = "Invalid cube in list file [" + ui.GetFileName("FROMLIST") + "]";
    throw IException(e, IException::User, msg, _FILEINFO_);
  }

  // Setup to propagate from the first input cube
  ProcessByLine p2;
  CubeAttributeInput inatt;

  int index = 0;
  if(ui.WasEntered("PROPLAB")) {
    bool match = false;
    QString fname = ui.GetFileName("PROPLAB");
    for(int i = 0; i < cubeList.size(); i++) {
      if(fname == cubeList[i].toString()) {
        index = i;
        match = true;
        break;
      }
    }
    if(!match) {
      QString msg = "FileName [" + ui.GetFileName("PROPLAB") +
                    "] to propagate labels from is not in the list file [" +
                    ui.GetFileName("FROMLIST") + "]";
      throw IException(IException::User, msg, _FILEINFO_);
    }
  }
  p2.SetInputCube(cubeList[index].toString(), inatt);


  // Create the output cube
  Cube *ocube = p2.SetOutputCube("TO", nsamps, nlines, nbands);
  p2.ClearInputCubes();

  p2.Progress()->SetText("Allocating cube");
  p2.StartProcess(NullBand);

  // Add the band bin group if necessary
  if(outBandBin.keywords() > 0) {
    ocube->putGroup(outBandBin);
  }
  p2.EndProcess();

  // Now loop and mosaic in each cube
  int sband = 1;
  for(int i = 0; i < cubeList.size(); i++) {
    ProcessMosaic m;
    m.SetBandBinMatch(false);
//.........这里部分代码省略.........

//.........这里部分代码省略.........
  cout << "nearest neighbor: " << chip.GetReadInterpolator() << endl;
  chip.SetReadInterpolator(Isis::Interpolator::BiLinearType);
  cout << "bilinear: " << chip.GetReadInterpolator() << endl;
  chip.SetReadInterpolator(Isis::Interpolator::CubicConvolutionType);
  cout << "cubic convolution: " << chip.GetReadInterpolator() << endl;

  cout << endl;
  cout << endl;
  cout << "Generate Errors:" << endl;
  Cube junkCube2;
  junkCube2.open("$base/testData/f319b18_ideal.cub");
  // 4 by 4 chip at samle 1000 line 500
  matchChip.TackCube(1, 1);
  matchChip.Load(junkCube2);

  cout << "Try to set interpolator to type 0 (Interpolator::None):" << endl;
  try {
    chip.SetReadInterpolator(Isis::Interpolator::None);
  }
  catch(IException &e) {
    ReportError(e.toString());
  }
  cout << "Try to set interpolator to type 3 (enum value not assigned):" << endl;
  try {
    chip.SetReadInterpolator((Isis::Interpolator::interpType) 3);
  }
  catch(IException &e) {
    ReportError(e.toString());
  }
  cout << "Try to set chip size with input parameter equal to 0:" << endl;
  try {
    newChip.SetSize(0, 1);
  }
  catch(IException &e) {
    ReportError(e.toString());
  }
  cout << "Try to load a cube that is not camera or map projection:" << endl;
  try {
    newChip.Load(junk, matchChip, junkCube);
  }
  catch(IException &e) {
    ReportError(e.toString());
  }
  cout << "Try to load a cube with a match cube that is not camera or map projection:" << endl;
  try {
    newChip.Load(junkCube, matchChip, junk);
  }
  catch(IException &e) {
    ReportError(e.toString());
  }
  cout << "Try to load a cube with match chip and cube that can not find at least 3 points for Affine Transformation:" << endl;
  try {
    newChip.Load(junkCube, matchChip, junkCube2);
  }
  catch(IException &e) {
    ReportError(e.toString());
  }
  cout << "Try to set valid range with larger number passed in as first parameter:" << endl;
  try {
    newChip.SetValidRange(4, 3);
  }
  catch(IException &e) {
    ReportError(e.toString());
  }
  cout << "Try to extract a sub-chip with samples or lines greater than original chip:" << endl;
  try {
    newChip.Extract(2, 5, 1, 1);
  }
  catch(IException &e) {
    ReportError(e.toString());
  }


  junk.close(true);// the "true" flag removes junk.cub from the /tmp/ directory
  junkCube.close(); // these cubes are kept in test data area, do not delete
  junkCube2.close();


#if 0
  try {
    junk.Open("/work2/janderso/moc/ab102401.lev1.cub");
    chip.TackCube(453.0, 567.0);
    chip.Load(junk);

    Cube junk2;
    junk2.Open("/work2/janderso/moc/ab102402.lev0.cub");
    Chip chip2(75, 70);
    chip2.TackCube(166.0, 567.0);
    chip2.Load(junk2, chip);

    chip.Write("junk3.cub");
    chip2.Write("junk4.cub");
  }
  catch(IException &e) {
    e.print();
  }
#endif

  return 0;
}

//.........这里部分代码省略.........
                     toString(readBytes) + "]" ;
        throw IException(IException::Io, msg, _FILEINFO_);
      }

      // Swap the bytes if necessary and convert any out of bounds pixels
      // to special pixels
      for(int samp = 0; samp < ns; samp++) {
        switch(pixelType) {
          case Isis::UnsignedByte:
            out[samp] = (double)((unsigned char *)in)[samp];
            break;
          case Isis::UnsignedWord:
            out[samp] = (double)swapper.UnsignedShortInt((unsigned short int *)in + samp);
            break;
          case Isis::SignedWord:
            out[samp] = (double)swapper.ShortInt((short int *)in + samp);
            break;
          case Isis::Real:
            out[samp] = (double)swapper.Float((float *)in + samp);
            break;
          default:
            break;
        }

        // Sets out to isis special pixel or leaves it if valid
        out[samp] = TestPixel(out[samp]);

        if(Isis::IsValidPixel(out[samp])) {
          out[samp] = mult * out[samp] + base;
        }
      } // End sample loop

      //Set the buffer position and write the line to the output file
      out.SetBasePosition(1, line + 1, band + 1);
      outCube.write(out);

      if(toInt(suffixItems[0]) != 0) {
        pos = fin.tellg();
        char *sideplaneData = new char[4*toInt(suffixItems[0])];
        fin.read(sideplaneData, 4 * toInt(suffixItems[0]));
        int suffixData = (int)swapper.Int((int *)sideplaneData);
        record[0] = line + 1;
        record[1] = band + 1;
        record[2] = suffixData;

        if(band < 96) {
          sideplaneVisTable += record;

          // set HIS pixels appropriately
          for(int samp = 0; samp < ns; samp++) {
            if(out[samp] >= 4095) {
              out[samp] = Isis::His;
            }
          }
        }
        else {
          record[1] = (band + 1) - 96;
          sideplaneIrTable += record;

          // set HIS pixels appropriately
          for(int samp = 0; samp < ns; samp++) {
            if(out[samp] + suffixData >= 4095) {
              out[samp] = Isis::His;
            }
          }
        }

        delete [] sideplaneData;

        // Check the last io
        if(!fin.good()) {
          QString msg = "Cannot read file [" + inFileName + "]. Position [" +
                       toString((int)pos) + "]. Byte count [" +
                       toString(4) + "]" ;
          throw IException(IException::Io, msg, _FILEINFO_);
        }
      }
    } // End band loop

    int backplaneSize = toInt(suffixItems[1]) * (4 * (ns + toInt(suffixItems[0])));
    fin.seekg(backplaneSize, ios_base::cur);

    // Check the last io
    if(!fin.good()) {
      QString msg = "Cannot read file [" + inFileName + "]. Position [" +
                   toString((int)pos) + "]. Byte count [" +
                   toString(4 * (4 * ns + 4)) + "]" ;
      throw IException(IException::Io, msg, _FILEINFO_);
    }

  } // End line loop

  outCube.write(sideplaneVisTable);
  outCube.write(sideplaneIrTable);

  // Close the file and clean up
  fin.close();
  outCube.close();
  delete [] in;
}

//.........这里部分代码省略.........
                    int l2 = ldim - 1;
                    for(int s = 0; s < sdim; s++) {
                        array[s][l1] = array[s][l2];
                    }
                }
                // Left Edge Reseau
                else if(toInt(type[res]) == 4) {
                    int s1 = 0;
                    int s2 = sdim - 1;
                    for(int l = 0; l < ldim; l++) {
                        array[s1][l] = array[s2][l];
                    }
                }
                // Right Edge Reseau
                else if(toInt(type[res]) == 6) {
                    int s1 = 0;
                    int s2 = sdim - 1;
                    for(int l = 0; l < ldim; l++) {
                        array[s2][l] = array[s1][l];
                    }
                }
                // Bottom Edge Reseau
                else if(toInt(type[res]) == 8) {
                    int l1 = 0;
                    int l2 = ldim - 1;
                    for(int s = 0; s < sdim; s++) {
                        array[s][l2] = array[s][l1];
                    }
                }
                // Walk top edge & replace data outside of 2devs with the avg
                for(int s = 0; s < sdim; s++) {
                    int l = 0;
                    double diff = fabs(array[s][l] - avg);
                    if(diff > (2 * sdev)) array[s][l] = avg;
                }
                // Walk bottom edge & replace data outside of 2devs with the avg
                for(int s = 0; s < sdim; s++) {
                    int l = ldim - 1;
                    double diff = fabs(array[s][l] - avg);
                    if(diff > (2 * sdev)) array[s][l] = avg;
                }
                // Walk left edge & replace data outside of 2devs with the avg
                for(int l = 0; l < ldim; l++) {
                    int s = 0;
                    double diff = fabs(array[s][l] - avg);
                    if(diff > (2 * sdev)) array[s][l] = avg;
                }
                // Walk right edge & replace data outside of 2devs with the avg
                for(int l = 0; l < ldim; l++) {
                    int s = sdim - 1;
                    double diff = fabs(array[s][l] - avg);
                    if(diff > (2 * sdev)) array[s][l] = avg;
                }
                srand(0);
                double dn, gdn1, gdn2;
                for(int l = 0; l < ldim; l++) {
                    int c = l * sdim;  //count
                    // Top Edge Reseau
                    if(toInt(type[res]) == 2 && l < (ldim / 2)) continue;
                    // Bottom Edge Reseau
                    if(toInt(type[res]) == 8 && l > (ldim / 2 + 1)) continue;
                    for(int s = 0; s < sdim; s++, c++) {
                        // Left Edge Reseau
                        if(toInt(type[res]) == 4 && s < (sdim / 2)) continue;
                        // Right Edge Reseau
                        if(toInt(type[res]) == 6 && s > (sdim / 2 + 1)) continue;
                        double sum = 0.0;
                        int gline1 = 0;
                        int gline2 = ldim - 1;
                        gdn1 = array[s][gline1];
                        gdn2 = array[s][gline2];

                        // Linear Interpolation to get pixel value
                        dn = gdn2 + (l - gline2) * (gdn1 - gdn2) / (gline1 - gline2);
                        sum += dn;

                        int gsamp1 = 0;
                        int gsamp2 = sdim - 1;
                        gdn1 = array[gsamp1][l];
                        gdn2 = array[gsamp2][l];

                        // Linear Interpolation to get pixel value
                        dn = gdn2 + (s - gsamp2) * (gdn1 - gdn2) / (gsamp1 - gsamp2);
                        sum += dn;
                        dn = sum / 2;
                        int rdm = rand();
                        double drandom = rdm / (double)RAND_MAX;
                        double offset = 0.0;
                        if(drandom < .333) offset = -1.0;
                        if(drandom > .666) offset = 1.0;
                        brick[c] = dn + offset;
                    }
                }
            }
        }
        cube.write(brick);
    }
    cube.close();

}

//.........这里部分代码省略.........
  c->SetUniversalGround(lat.degrees(), lon.degrees(), radius);
  c->SetGround(SurfacePoint(lat, lon, Distance(radius, Distance::Meters)));
  cout << "Has intersection " << c->HasSurfaceIntersection() << endl;
  cout << "Latitude = " << c->UniversalLatitude() << endl;
  cout << "Longitude = " << c->UniversalLongitude() << endl;
  cout << "Radius = " << c->LocalRadius().meters() << endl;
  double p[3];
  c->Coordinate(p);
  cout << "Point = " << setprecision(4)
       << p[0] << " " << p[1] << " " << p[2] << endl << endl;

  cout << "Test Forward/Reverse Camera Calculations At Center Of Image..."
       << endl;
  sample = c->Samples() / 2.0;
  line = c->Lines() / 2.0;
  cout << "Sample = " << setprecision(3) << sample << endl;
  cout << "Line = " << line << endl;
  cout << "SetImage (sample, line): " << c->SetImage(sample, line) << endl;
  cout << "Latitude = " << c->UniversalLatitude() << endl;
  cout << "Longitude = " << c->UniversalLongitude() << endl;
  cout << "Radius = " << c->LocalRadius().meters() << endl;
  c->Coordinate(p);
  cout << "Point = " << setprecision(4)
       << p[0] << " " << p[1] << " " << p[2] << endl;
  cout << "SetUniversalGround (lat, lon, radius): "
       << c->SetUniversalGround(c->UniversalLatitude(), c->UniversalLongitude(),
                                c->LocalRadius().meters())
       << endl;
  cout << "Sample = " << c->Sample() << endl;
  cout << "Line = " << c->Line() << endl << endl;

  cout << endl << "/---------- Test Polar Boundary Conditions" << endl;
  inputFile = "$clementine1/testData/lub5992r.292.lev1.phot.cub";
  cube.close();
  cube.open(inputFile);
  pvl = *cube.label();
  Camera *cam2 = CameraFactory::Create(cube);
  cube.close();

  cout << endl;
  cout << "Camera* from: " << inputFile << endl;
  ifovOffsets = cam2->PixelIfovOffsets();
  cout << "Pixel Ifov: " << endl;
  foreach (QPointF offset, ifovOffsets) {
    cout << offset.x() << " , " << offset.y() << endl;
  }
  cout << "Basic Mapping: " << endl;
  Pvl camMap;
  cam2->BasicMapping(camMap);

  double minLat = camMap.findGroup("Mapping")["MinimumLatitude"];
  minLat *= 100;
  minLat = round(minLat);
  minLat /= 100;
  camMap.findGroup("Mapping")["MinimumLatitude"] = toString(minLat);

  double pixRes = camMap.findGroup("Mapping")["PixelResolution"];
  pixRes *= 100;
  pixRes = round(pixRes);
  pixRes /= 100;
  camMap.findGroup("Mapping")["PixelResolution"] = toString(pixRes);

  double minLon = camMap.findGroup("Mapping")["MinimumLongitude"];
  minLon *= 100000000000.0;
  minLon = round(minLon);
  minLon /= 100000000000.0;

void IsisMain() {

  Preference::Preferences(true);

  cout << "Testing Isis::ProcessMapMosaic Class ... " << endl;

  // Create the temp parent cube
  FileList cubes;
  cubes.read("unitTest.lis");

  cout << "Testing Mosaic 1" << endl;
  ProcessMapMosaic m1;
  CubeAttributeOutput oAtt;
  ProcessMosaic::ImageOverlay priority = ProcessMapMosaic::PlaceImagesOnTop;
  m1.SetBandBinMatch(false);
  m1.SetOutputCube(cubes, oAtt, "./unitTest.cub");

  //set priority
  m1.SetImageOverlay(priority);

  for(int i = 0; i < cubes.size(); i++) {
    if(m1.StartProcess(cubes[i].toString())) {
      cout << cubes[i].toString() << " is inside the mosaic" << endl;
    }
    else {
      cout << cubes[i].toString() << " is outside the mosaic" << endl;
    }
  }

  m1.EndProcess();
  cout << "Mosaic label: " << endl;

  Pvl labels("./unitTest.cub");
  cout << labels << endl;

  remove("./unitTest.cub");

  cout << "Testing Mosaic 2" << endl;
  ProcessMapMosaic m2;
  m2.SetBandBinMatch(false);
  m2.SetOutputCube(cubes, -6, -4, 29, 31, oAtt, "./unitTest.cub");

  //set priority
  m2.SetImageOverlay(priority);

  for(int i = 0; i < cubes.size(); i++) {
    if(m2.StartProcess(cubes[i].toString())) {
      cout << cubes[i].toString() << " is inside the mosaic" << endl;
    }
    else {
      cout << cubes[i].toString() << " is outside the mosaic" << endl;
    }
  }

  m2.EndProcess();
  cout << "Mosaic label: " << endl;

  labels.clear();
  labels.read("./unitTest.cub");
  cout << labels << endl;

  Cube tmp;
  tmp.open("./unitTest.cub");
  LineManager lm(tmp);
  lm.SetLine(1, 1);

  while(!lm.end()) {
    tmp.read(lm);
    cout << "Mosaic Data: " << lm[lm.SampleDimension()/4] << '\t' <<
              lm[lm.SampleDimension()/2] << '\t' <<
              lm[(3*lm.SampleDimension())/4] << endl;
    lm++;
  }

  tmp.close();
  remove("./unitTest.cub");  // Create the temp parent cube

  cout << endl << "Testing Mosaic where the input (x, y) is negative,"
          " according to the output cube." << endl;
  QString inputFile = "./unitTest1.cub";
  Cube inCube;
  inCube.open(inputFile);
  PvlGroup mapGroup = inCube.label()->findGroup("Mapping", Pvl::Traverse);

  mapGroup.addKeyword(PvlKeyword("MinimumLatitude",  toString(-4.9)), Pvl::Replace);
  mapGroup.addKeyword(PvlKeyword("MaximumLatitude",  toString(-4.7)), Pvl::Replace);
  mapGroup.addKeyword(PvlKeyword("MinimumLongitude", toString(30.7)), Pvl::Replace);
  mapGroup.addKeyword(PvlKeyword("MaximumLongitude", toString(31)), Pvl::Replace);
  
  inCube.close();
  CubeAttributeOutput oAtt2( FileName("./unitTest3.cub") );
  ProcessMapMosaic m3;
  
  m3.SetBandBinMatch(false);
  m3.SetOutputCube(inputFile, mapGroup, oAtt2, "./unitTest3.cub");

  //set priority
  m3.SetImageOverlay(priority);
  m3.SetHighSaturationFlag(false);
  m3.SetLowSaturationFlag(false);
//.........这里部分代码省略.........

//.........这里部分代码省略.........

  // If a TO parameter was specified, create DEM with errors
  if (ui.WasEntered("TO")) {
    //  Create the output DEM
    cout << "\nCreating output DEM from " << bmf.size() << " points.\n";
    Process  p;
    Cube *icube = p.SetInputCube("FROM");
    Cube *ocube = p.SetOutputCube("TO", icube->sampleCount(),
                                  icube->lineCount(), 3);
    p.ClearInputCubes();

    int boxsize = ui.GetInteger("BOXSIZE");
    double plotdist = ui.GetDouble("PLOTDIST");

    TileManager dem(*ocube), eigen(*ocube), stErr(*ocube);
    dem.SetTile(1, 1);      //  DEM Data/elevation
    stErr.SetTile(1, 2);    //  Error in stereo computation
    eigen.SetTile(1, 3);    //  Eigenvalue of the solution

    int nBTiles(eigen.Tiles()/3);  // Total tiles / 3 bands

    prog.SetText("Creating DEM");
    prog.SetMaximumSteps(nBTiles);
    prog.CheckStatus();

    Statistics stAng;
    while ( !eigen.end() ) {   // Must use the last band for this!!
      PointPlot tm = for_each(bmf.begin(), bmf.end(), PointPlot(dem, plotdist));
      tm.FillPoints(*lhCamera, *rhCamera, boxsize, dem, stErr, eigen, &stAng);

      ocube->write(dem);
      ocube->write(stErr);
      ocube->write(eigen);

      dem.next();
      stErr.next();
      eigen.next();

      prog.CheckStatus();
    }

    //  Report Stereo separation angles
    PvlGroup stresultsPvl("StereoSeparationAngle");
    stresultsPvl += PvlKeyword("Minimum", toString(stAng.Minimum()), "deg");
    stresultsPvl += PvlKeyword("Average", toString(stAng.Average()), "deg");
    stresultsPvl += PvlKeyword("Maximum", toString(stAng.Maximum()), "deg");
    stresultsPvl += PvlKeyword("StandardDeviation", toString(stAng.StandardDeviation()), "deg");
    Application::Log(stresultsPvl);

    // Update the label with BandBin keywords
    PvlKeyword filter("FilterName", "Elevation", "meters");
    filter.addValue("ElevationError", "meters");
    filter.addValue("GoodnessOfFit", "unitless");
    PvlKeyword center("Center", "1.0");
    center.addValue("1.0");
    center.addValue("1.0");

    PvlGroup &bandbin = ocube->label()->findGroup("BandBin", PvlObject::Traverse);
    bandbin.addKeyword(filter, PvlContainer::Replace);
    bandbin.addKeyword(center, PvlContainer::Replace);
    center.setName("Width");
    bandbin.addKeyword(center, PvlContainer::Replace);


    p.EndProcess();
  }

  // If a cnet file was entered, write the ControlNet pvl to the file
  if (ui.WasEntered("ONET")) {
    WriteCnet(ui.GetFileName("ONET"), bmf, lhCamera->target()->name(), serialLeft,
              serialRight);
  }

  // Create output data
  PvlGroup totalPointsPvl("Totals");
  totalPointsPvl += PvlKeyword("AttemptedPoints", toString(numAttemptedInitialPoints));
  totalPointsPvl += PvlKeyword("InitialSuccesses", toString(numOrigPoints));
  totalPointsPvl += PvlKeyword("GrowSuccesses", toString(passpix2));
  totalPointsPvl += PvlKeyword("ResultingPoints", toString(bmf.size()));

  Application::Log(totalPointsPvl);

  Pvl arPvl = matcher.RegistrationStatistics();
  PvlGroup smtkresultsPvl("SmtkResults");
  smtkresultsPvl += PvlKeyword("SpiceOffImage", toString(matcher.OffImageErrorCount()));
  smtkresultsPvl += PvlKeyword("SpiceDistanceError", toString(matcher.SpiceErrorCount()));
  arPvl.addGroup(smtkresultsPvl);

  for(int i = 0; i < arPvl.groups(); i++) {
    Application::Log(arPvl.group(i));
  }

  // add the auto registration information to print.prt
  PvlGroup autoRegTemplate = matcher.RegTemplate();
  Application::Log(autoRegTemplate);

  // Don't need the cubes opened anymore
  lhImage.close();
  rhImage.close();
}

//.........这里部分代码省略.........

  //record first fiducial measurement in the table
  recordFid[0] = 0;
  recordFid[1] = sampleInitial;
  recordFid[2] = lineInitial;
  tableFid += recordFid;
  for (s= sampleInitial, l=lineInitial, fidn=0;  s<panS;  s+=averageSamples, fidn++) {
     //corrections for half spacing of center fiducials
     if (fidn == 22) s -= averageSamples/2.0;
     if (fidn == 23) s -= averageSamples/2.0;

     //look for the bottom fiducial
     searchS.TackCube(s,l+averageLines);
     searchS.Load(panCube, 0, scale);
     *arS->SearchChip()   = searchS;
     regStatus = arS->Register();
     if (regStatus == AutoReg::SuccessPixel) {
       inputChip.TackCube(arS->CubeSample(), arS->CubeLine());
       inputChip.Load(panCube,0,1);
       inputChip.SetCubePosition(arS->CubeSample(), arS->CubeLine());
     }
     else {  //if autoreg is unsuccessful, a larger window will be used
       inputChip.TackCube(s, l+averageLines);
       inputChip.Load(panCube, 0, 1);
       inputChip.SetCubePosition(s, l+averageLines);
     }
     centroid.selectAdaptive(&inputChip, &selectionChip);  //continuous dynamic range selection
     //elliptical trimming/smoothing... if this fails move on
     if (centroid.elipticalReduction(&selectionChip, 95, play, 2000) != 0 ) {      
       //center of mass to reduce selection to a single measure
       centroid.centerOfMass(&selectionChip, &sample, &line);      
       inputChip.SetChipPosition(sample, line);
       sample = inputChip.CubeSample();
       line   = inputChip.CubeLine();
       recordFid[0] = fidn*2+1;
       recordFid[1] = sample;
       recordFid[2] = line;
       tableFid += recordFid;
     }
     progress.CheckStatus();

     //look for the top fiducial
     if (s == sampleInitial) //first time through the loop?
       continue;  //then the top fiducial was already found
     searchS.TackCube(s, l);
     searchS.Load(panCube, 0, scale);
     *arS->SearchChip()   = searchS;
     regStatus = arS->Register();
     if (regStatus == AutoReg::SuccessPixel) {
       inputChip.TackCube(arS->CubeSample(), arS->CubeLine());
       inputChip.Load(panCube, 0, 1);
       inputChip.SetCubePosition(arS->CubeSample(), arS->CubeLine());
     }
     else {  //if autoreg is unsuccessful, a larger window will be used
       inputChip.TackCube(s, l);
       inputChip.Load(panCube, 0, 1);
       inputChip.SetCubePosition(s, l);
     }
     centroid.selectAdaptive(&inputChip, &selectionChip);//continuous dynamic range selection
     //inputChip.Write("inputTemp.cub");//debug
     //selectionChip.Write("selectionTemp.cub");//debug
     //elliptical trimming/smoothing... if this fails move on
     if (centroid.elipticalReduction(&selectionChip, 95, play, 2000) !=0) {    
       //center of mass to reduce selection to a single measure
       centroid.centerOfMass(&selectionChip, &sample, &line);  
       inputChip.SetChipPosition(sample, line);
       //when finding the top fiducial both s and l are refined for a successful measurement, 
       //  this will help follow trends in the scaned image
       s = inputChip.CubeSample(); 
       l = inputChip.CubeLine();
       recordFid[0] = fidn*2;
       recordFid[1] = s;
       recordFid[2] = l;
       tableFid += recordFid;
     }
     progress.CheckStatus();
  }

  panCube.write(tableFid);
  //close the new cube
  panCube.close(false);
  panCube.open(ui.GetFileName("FROM"),"rw");
 
  delete spPos;
  delete spRot;

  //now instantiate a camera to make sure all of this is working
  ApolloPanoramicCamera* cam = (ApolloPanoramicCamera*)(panCube.camera());
  //log the residual report from interior orientation 
  PvlGroup residualStats("InteriorOrientationStats");
  residualStats += PvlKeyword("FiducialsFound",  toString(tableFid.Records()));
  residualStats += PvlKeyword("ResidualMax",  toString(cam->intOriResidualMax()),"pixels");
  residualStats += PvlKeyword("ResidualMean", toString(cam->intOriResidualMean()),"pixels");
  residualStats += PvlKeyword("ResidualStdev", toString(cam->intOriResidualStdev()),"pixels");

  Application::Log( residualStats ); 


  return;
}

//.........这里部分代码省略.........
                       (poly.validLineDim() * 2) - 3.0) /
                       ui.GetInteger("NUMVERTICES")));
    if (sinc < 1.0 || linc < 1.0)
      sinc = linc = 1.0;
  }
  else if (incType.UpCase() == "LINCSINC"){
    sinc = ui.GetInteger("SINC");
    linc = ui.GetInteger("LINC");
  }
  else {
    string msg = "Invalid INCTYPE option[" + incType + "]";
    throw IException(IException::Programmer, msg, _FILEINFO_);
  }

  bool precision = ui.GetBoolean("INCREASEPRECISION");
  try {
    poly.Create(cube, sinc, linc, 1, 1, 0, 0, 1, precision);
  }
  catch (IException &e) {
    QString msg = "Cannot generate polygon for [" + ui.GetFileName("FROM") + "]";
    throw IException(e, IException::User, msg, _FILEINFO_);
  }

  if(ui.GetBoolean("TESTXY")) {
    Pvl cubeLab(ui.GetFileName("FROM"));
    PvlGroup inst = cubeLab.findGroup("Instrument", Pvl::Traverse);
    QString target = inst["TargetName"];
    PvlGroup radii = Projection::TargetRadii(target);

    Pvl map(ui.GetFileName("MAP"));
    PvlGroup &mapping = map.findGroup("MAPPING");

    if(!mapping.hasKeyword("TargetName"))
      mapping += Isis::PvlKeyword("TargetName", target);
    if(!mapping.hasKeyword("EquatorialRadius"))
      mapping += Isis::PvlKeyword("EquatorialRadius", (QString)radii["EquatorialRadius"]);
    if(!mapping.hasKeyword("PolarRadius"))
      mapping += Isis::PvlKeyword("PolarRadius", (QString)radii["PolarRadius"]);
    if(!mapping.hasKeyword("LatitudeType"))
      mapping += Isis::PvlKeyword("LatitudeType", "Planetocentric");
    if(!mapping.hasKeyword("LongitudeDirection"))
      mapping += Isis::PvlKeyword("LongitudeDirection", "PositiveEast");
    if(!mapping.hasKeyword("LongitudeDomain"))
      mapping += Isis::PvlKeyword("LongitudeDomain", "360");
    if(!mapping.hasKeyword("CenterLatitude"))
      mapping += Isis::PvlKeyword("CenterLatitude", "0");
    if(!mapping.hasKeyword("CenterLongitude"))
      mapping += Isis::PvlKeyword("CenterLongitude", "0");

    sinc = poly.getSinc();
    linc = poly.getLinc();
    bool polygonGenerated = false;
    while (!polygonGenerated) {
      Projection *proj = NULL;
      geos::geom::MultiPolygon *xyPoly = NULL;

      try {
        proj = ProjectionFactory::Create(map, true);
        xyPoly = PolygonTools::LatLonToXY(*poly.Polys(), proj);

        polygonGenerated = true;
      }
      catch (IException &e) {
        if (precision && sinc > 1 && linc > 1) {
          sinc = sinc * 2 / 3;
          linc = linc * 2 / 3;
          poly.Create(cube, sinc, linc);
        }
        else {
          delete proj;
          delete xyPoly;
          e.print(); // This should be a NAIF error
          QString msg = "Cannot calculate XY for [";
          msg += ui.GetFileName("FROM") + "]";
          throw IException(e, IException::User, msg, _FILEINFO_);
        }
      }

      delete proj;
      delete xyPoly;
    }
  }

  cube.deleteBlob("Polygon", sn);
  cube.write(poly);

  if(precision) {
    PvlGroup results("Results");
    results.addKeyword(PvlKeyword("SINC", toString(sinc)));
    results.addKeyword(PvlKeyword("LINC", toString(linc)));
    Application::Log(results);
  }

  Process p;
  p.SetInputCube("FROM");
  p.WriteHistory(cube);

  cube.close();
  prog.CheckStatus();
}

  /**
   * Set the output cube to specified file name and specified input images
   * and output attributes
   */
  Isis::Cube *ProcessMapMosaic::SetOutputCube(FileList &propagationCubes, CubeAttributeOutput &oAtt,
      const QString &mosaicFile) {
    int bands = 0;
    double xmin = DBL_MAX;
    double xmax = -DBL_MAX;
    double ymin = DBL_MAX;
    double ymax = -DBL_MAX;
    double slat = DBL_MAX;
    double elat = -DBL_MAX;
    double slon = DBL_MAX;
    double elon = -DBL_MAX;

    Projection *proj = NULL;

    if (propagationCubes.size() < 1) {
      QString msg = "The list does not contain any data";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }

    for (int i = 0; i < propagationCubes.size(); i++) {
      // Open the cube and get the maximum number of band in all cubes
      Cube cube;
      cube.open(propagationCubes[i].toString());
      bands = max(bands, cube.bandCount());

      // See if the cube has a projection and make sure it matches
      // previous input cubes
      Projection *projNew =
          Isis::ProjectionFactory::CreateFromCube(*(cube.label()));
      if ((proj != NULL) && (*proj != *projNew)) {
        QString msg = "Mapping groups do not match between cubes [" +
                     propagationCubes[0].toString() + "] and [" + propagationCubes[i].toString() + "]";
        throw IException(IException::User, msg, _FILEINFO_);
      }

      // Figure out the x/y range as it may be needed later
      double x = projNew->ToProjectionX(0.5);
      double y = projNew->ToProjectionY(0.5);
      if (x < xmin) xmin = x;
      if (y < ymin) ymin = y;
      if (x > xmax) xmax = x;
      if (y > ymax) ymax = y;

      x = projNew->ToProjectionX(cube.sampleCount() + 0.5);
      y = projNew->ToProjectionY(cube.lineCount() + 0.5);
      if (x < xmin) xmin = x;
      if (y < ymin) ymin = y;
      if (x > xmax) xmax = x;
      if (y > ymax) ymax = y;

      slat = min(slat, projNew->MinimumLatitude());
      elat = max(elat, projNew->MaximumLatitude());
      slon = min(slon, projNew->MinimumLongitude());
      elon = max(elon, projNew->MaximumLongitude());

      // Cleanup
      cube.close();
      if (proj) delete proj;
      proj = projNew;
    }

    if (proj) delete proj;

    return SetOutputCube(propagationCubes[0].toString(), xmin, xmax, ymin, ymax,
                         slat, elat, slon, elon, bands, oAtt, mosaicFile);
  }

//.........这里部分代码省略.........
  const int nlBefore = 2;
  const int nlAfter = 2;
  const int nlavg = 4;
  const double fixFactor = 1.5;

  // Loop for each line to fix
  for(unsigned int i = 0; i < fixList.size(); i++) {
    // We will examine a window of lines around the wago change
    int centerLine = fixList[i];
    int sl = centerLine - nlBefore - nlavg;
    int el = centerLine + nlAfter + nlavg;

    // Don't do anything if the window is outside the image
    if(sl < 1) continue;
    if(el < 1) continue;
    if(sl > nl) continue;
    if(el > nl) continue;

    // Find the closest non-zero output line average before the wago line
    int slBefore = sl;
    int elBefore = sl + nlavg;
    int indexBefore = -1;
    for(int line = elBefore; line >= slBefore; line--) {
      if(gbl::outLineAvg[line-1] != 0.0) {
        indexBefore = line - 1;
        break;
      }
    }
    if(indexBefore < 0) continue;
    double oavgBefore = gbl::outLineAvg[indexBefore];

    // Find the closest non-zero output line average before the wago line
    int slAfter = el - nlavg;
    int elAfter = el;
    int indexAfter = -1;
    for(int line = slAfter; line <= elAfter; line++) {
      if(gbl::outLineAvg[line-1] != 0.0) {
        indexAfter = line - 1;
        break;
      }
    }
    if(indexAfter < 0) continue;
    double oavgAfter = gbl::outLineAvg[indexAfter];

    // Get the corresponding input averages and compute net WAGO change
    // Don't do anything if the net change is invalid
    double iavgBefore = gbl::inLineAvg[indexBefore];
    double iavgAfter  = gbl::inLineAvg[indexAfter];

    double sum = (iavgBefore / iavgAfter) / (oavgBefore / oavgAfter);
    if(abs(1.0 - sum) < 0.05) continue;

    // Prep to fix the lines
    sl = centerLine - nlBefore;
    el = centerLine + nlAfter;
    int nlFix = el - sl + 1;
    double fixinc = (oavgAfter - oavgBefore) / (nlFix + 1);
    double base = oavgBefore;
    double avgTol = abs(fixinc * fixFactor);

    // Loop and fix each one
    for(int line = sl; line <= el; line++) {
      base += fixinc;
      double avg = base - gbl::outLineAvg[line-1];

      // Do we need to fix this line?
      if(abs(avg) <= avgTol) continue;

      // Read the line
      lbuf.SetLine(line);
      fix.read(lbuf);

      // Null it
      if(gbl::nullWago) {
        for(int samp = 0 ; samp < lbuf.size(); samp++) {
          lbuf[samp] = Null;
        }
        outLineAvg[line-1] = 0.0;
      }

      // or repair it
      else {
        avg = outLineAvg[line-1];
        outLineAvg[line-1] = base;

        for(int samp = 0; samp < lbuf.size(); samp++) {
          if(IsValidPixel(lbuf[samp])) {
            lbuf[samp] = lbuf[samp] / avg * base;
          }
        }
      }

      // Write the line
      fix.write(lbuf);
    }
  }

  // Cleanup
  fix.close();
}