C++ To::bscrunch方法代码示例

void Pulsar::SquareWave::get_transitions (const Profile* profile,
					  vector<unsigned>& up,
					  vector<unsigned>& down)
{
  unsigned nbin = profile->get_nbin();
  unsigned scrunch = 1;

  Reference::To<Profile> clone;
  if (use_nbin && nbin > use_nbin) {
    clone = profile->clone();
    scrunch = nbin/use_nbin;
    clone->bscrunch(scrunch);
    profile = clone;
    nbin = profile->get_nbin();
  }

  unsigned offset = (unsigned) (risetime * nbin);

cerr << "nbin=" << nbin << " offset=" << offset << endl;

  // differentiate the profile
  Reference::To<Profile> difference = differentiate (profile, offset);

  float* amps = difference->get_amps();

  // find the phase window in which the mean is closest to zero
  BaselineWindow window;
  window.set_find_mean (0.0);
  window.get_smooth()->set_turns (0.2);
  float zero = window.find_phase (nbin, amps);

  // get the noise statistics of the zero mean region
  double mean = 0;
  double variance = 0;
  difference->stats (zero, &mean, &variance);

cerr << "mean=" << mean << " rms=" << sqrt(variance) << endl;

  // check that the mean is actually zero
  double rms = sqrt(variance);
  if (mean > rms)
    throw Error (InvalidState, "Pulsar::SquareWave::get_transitions",
		 "mean=%lf > rms=%lf", mean, rms);

  float cutoff = threshold * rms;

  find_transitions (nbin, amps, up, cutoff);
  find_transitions (nbin, amps, down, -cutoff);
}

void psrspa::matched_finder ( const Archive* arch )
{
  string name = arch->get_filename ();
  for ( unsigned isub = 0; isub < arch->get_nsubint(); isub++ )
  {
    Reference::To<Profile> profile = arch->get_Profile(isub,0,0)->clone();

    while (profile->get_nbin() > 256)
    {
      finder->set_Profile( profile );  // might finder optimize on &profile?

      PhaseWeight weight;
      finder->get_weight( &weight );

      matched_report (name, isub, weight, *profile);

      profile->bscrunch (2);
    }
  }
}

//.........这里部分代码省略.........
	command += " -f ";
	command += optarg;
	break;
	
      case 'n':

	reflections.add_reflection( optarg[0] );

	command += " -n ";
	command += optarg;
	break;

      case 'o':
	new_cfreq = true;
	if (sscanf(optarg, "%lf", &new_fr) != 1) {
	  cout << "That is not a valid centre frequency" << endl;
	  return -1;
	}
	command += " -o ";
	command += optarg;
	break;
      case 't':
	tscr = true;
	if (sscanf(optarg, "%d", &tscr_fac) != 1) {
	  cout << "That is not a valid tscrunch factor" << endl;
	  return -1;
	}
	command += " -t ";
	command += optarg;
	break;
      case 'b':
	bscr = true;
	if (sscanf(optarg, "%d", &bscr_fac) != 1) {
	  cout << "That is not a valid bscrunch factor" << endl;
	  return -1;
	}
	if (bscr_fac <= 0) {
	  cout << "That is not a valid bscrunch factor" << endl;
	  return -1;
	}
	command += " -b ";
	command += optarg;
	break;
      case 'd':
	newdm = true;
	if (sscanf(optarg, "%lf", &dm) != 1) {
	  cout << "That is not a valid dispersion measure" << endl;
	  return -1;
	}
	command += " -d ";
	command += optarg;
	break;
      case 'D':
	dedisperse = true;
	command += " -D ";
	break;
      case 'R':
	if (sscanf(optarg, "%lf", &rm) != 1) {
	  cout << "That is not a valid rotation measure" << endl;
	  return -1;
	}
	newrm = true;
	defaraday = true;
	command += " -R ";
	command += optarg;
	break;

int main (int argc, char** argv) try
{
  // the multiple component model
  ComponentModel model;

  string model_filename_in;
  string model_filename_out = "paas.m";
  string details_filename   = "paas.txt";
  string std_filename       = "paas.std";

  bool fit = false;
  vector<string> new_components;
  string fit_flags;

  int bscrunch = 1;

  string pgdev;

  bool line_plot=false;
  bool interactive = false;
  bool centre_model = false, rotate_peak=false;
  float rotate_amount = 0.0;
  bool align = false;

  const char* args = "hb:r:w:c:fF:it:d:Dl:j:Ws:CpR:a";
  int c;

  while ((c = getopt(argc, argv, args)) != -1)
    switch (c) {

    case 'h':
      usage ();
      return 0;

    case 'r':
      model_filename_in = optarg;
      break;
      
    case 'w':
      model_filename_out = optarg;
      break;
      
    case 'f':
      fit = true;
      break;

    case 'W':
      model.set_fit_derivative (true);
      break;

    case 'c':
      new_components.push_back(optarg);
      break;

    case 'F':
      fit_flags = optarg;
      break;

    case 'b':
      bscrunch = atoi (optarg);
      break;

    case 't':
      model.set_threshold( atof(optarg) );
      break;

    case 'd':
      pgdev = optarg;
      break;

    case 'i':
      interactive = true;
    case 'D':
      pgdev = "/xs";
      break;
      
    case 'l':
      line_plot = true;
      break;

    case 's':
      std_filename = optarg;
      break;

    case 'C':
      centre_model = true;
      break;

    case 'p':
      rotate_peak = true;
      break;

    case 'R':
      rotate_amount = atof(optarg);
      break;

    case 'a':
      align = true;
      break;

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

void psrspa::create_histograms ( Reference::To<Archive> archive )
{
  if ( verbose )
    cerr << "psrspa::create_histograms entered" << endl;
  // ensure Stokes parameters if creating polarisation histograms
  if ( create_polar_degree || create_polar_angle )
  {
    archive->convert_state ( Signal::Stokes );
    if ( verbose )
      cerr << "psrspa::create_histograms converted state of the archive to Stokes" << endl;
  }

  // auxillary vectors
  //vector< Estimate<float > > aux_vec_f;
  vector< Estimate<double > > aux_vec_d;

  // Full Stokes profile
  Reference::To<PolnProfile> profile;
  // Polarized flux
  Reference::To<Profile> P;
  P = new Profile;
  // Total flux
  float *T_amps = new float [ nbin ];
  float *P_amps = new float [ nbin ];

  unsigned bin_min, bin_max;

  if ( verbose && max_bscrunch > 1 )
    cerr << "psrspa::create_histograms entering the bscrunch loop for the " << archive->get_filename () << " archive" << endl;
  // TODO Uh, this should be fixed up. The first idea was to enable the phase resolved histograms to be bscrunch-aware as well, but I think it doesn't make much sense so I decided to keep only the max flux bscrunch aware. This can be done much more neatly probably in such a case
  for ( current_bscrunch = 1 ; current_bscrunch <= max_bscrunch ; current_bscrunch *= 2 )
  {
    // in each passage, we bscrunch by a factor of 2
    if ( current_bscrunch > 1 )
    {
      if ( verbose )
	cerr << "psrspa::create_histograms bscrunching the archive " << archive->get_filename () << " by a factor of 2" << endl;
      archive->bscrunch ( 2 );
    }
    if ( verbose )
      cerr << "psrspa::create_histograms entering the loop through subints of the " << archive->get_filename () << " archive" << endl;
    // loop through subints
    for ( unsigned isub = 0; isub < archive->get_nsubint (); isub++ )
    {
      if ( verbose )
	cerr << "psrspa::create_histograms creating necessary profiles for subint " << isub << " of " << archive->get_filename () << endl;
      if ( create_polar_angle || create_polar_degree && current_bscrunch == 1 )
      {
	profile = archive->get_Integration(isub)->new_PolnProfile(0);
	if ( verbose )
	  cerr << "psrspa::create_histograms retrieved PolnProfile for subint " << isub << " of " << archive->get_filename () << endl;
      }
      if ( create_polar_angle && current_bscrunch == 1 )
      {
	profile->get_orientation ( aux_vec_d, 0 );
	if ( verbose )
	  cerr << "psrspa::create_histograms retrieved polarisation angle for subint " << isub << " of " << archive->get_filename () << endl;
      }
      if ( create_polar_degree || create_flux || find_max_amp_in_range )
      {
	stats.set_profile ( archive->get_Profile ( isub, 0, 0 ) );
	b_sigma = sqrt ( stats.get_baseline_variance ().get_value () );
	T_amps = archive->get_Profile ( isub, 0, 0 )->get_amps (); 
	if ( verbose )
	  cerr << "psrspa::create_histograms retrieved total flux amps for subint " << isub << " of " << archive->get_filename () << endl;
	if ( create_polar_degree && current_bscrunch == 1 )
	{
	  profile->get_polarized ( P );
	  if ( verbose )
	    cerr << "psrspa::create_histograms retrieved polarized flux profile for subint "  << isub << " of " << archive->get_filename () << endl;
	  P_amps = P->get_amps ();
	  if ( verbose )
	    cerr << "psrspa::create_histograms retrieved polarized flux amps for subint " << isub << " of " << archive->get_filename () << endl;
	}
      }

      if ( verbose )
	cerr << "psrspa::create_histograms looping through the provided phase ranges for subint " << isub << " of " << archive->get_filename () << endl;
      unsigned curr_hist = 0;
      // loop through phase ranges
      for ( unsigned irange = 0; irange < phase_range.size () ; irange++ )
      {
	if ( irange%2 == 0)
	{
	  bin_min = unsigned( floor ( phase_range[irange] * float(nbin / current_bscrunch ) + 0.5 ) );
	  if ( verbose ) 
	    cerr << "psrspa::create_histograms set minimal bin to " << bin_min << endl;
	}
	else
	{
	  bin_max = unsigned( floor ( phase_range[irange] * float(nbin / current_bscrunch ) + 0.5 ) );
	  if ( bin_max == nbin )
	    bin_max = nbin - 1 ;
	  if ( verbose ) 
	    cerr << "psrspa::create_histograms set maximum bin to " << bin_max << endl;
	  // loop through bins in the given phase range
	  for ( unsigned ibin = bin_min; ibin <= bin_max; ibin++ )
	  {
	    if ( create_polar_angle && current_bscrunch == 1 )
	    {
//.........这里部分代码省略.........