C++ Molecule::getProperty方法代码示例

void mergeDRFiles(vector<String>& names, string& output_file, Size& best_k, string& e_property, double& score_cutoff, double& score_cuton)
{

	DockResultFile* output = new DockResultFile(output_file, ios::out);
	bool sort_by_scores = 1;
	if (e_property == "") sort_by_scores = 0;

	vector<Result*> new_results;

	/// First of all, copy Result data
	map<Result::Method, Result*> result_map;
	for (Size file = 0; file < names.size(); file++)
	{
		DockResultFile* input = new DockResultFile(names[file]);

		const vector<Result*>* results = input->getResults();
		for (Size i = 0; i < results->size(); i++)
		{
			map<Result::Method, Result*>::iterator it = result_map.find((*results)[i]->getMethod());
			if (it == result_map.end())
			{
				Result* result_copy = new Result(*(*results)[i]);
				if (!sort_by_scores) output->addResult(result_copy);
				else new_results.push_back(result_copy);
				result_map.insert(make_pair(result_copy->getMethod(), result_copy));
			}
			else
			{
				*it->second += *(*results)[i];
			}
		}

		input->close();
		delete input;
	}


	if (e_property != "")
	{
		e_property = "score_"+new_results.back()->getMethodString();
	}

	/// If no sorting is desired, iterate over all input-files and write each input-molecules to output-file
	if (!sort_by_scores)
	{
		output->disableAutomaticResultCreation();

		for (Size file = 0; file < names.size(); file++)
		{
			GenericMolFile* input = MolFileFactory::open(names[file]);

			int mol_no = 0;
			for (Molecule* mol = input->read(); mol; mol = input->read(), mol_no++)
			{
				*output << *mol;
				delete mol;
				Log.level(20) << "\r" << names[file] << " : " << mol_no+1;
				Log.flush();
			}
			Log.level(20)<<endl;
			Log.flush();

			input->close();
			delete input;
		}
	}

	/// If sorting is desired, iterate over all input-files and save each input-molecules to a map.
	/// Then write all FlexibleMolecules in this map to the output file and adapt the Result objects.
	else
	{
		multimap < double, FlexibleMolecule* > compounds; // map containing score and conformation-ID
		set < String > IDs; // IDs of the base-conformations

		for (Size file = 0; file < names.size(); file++)
		{
			DockResultFile* input = new DockResultFile(names[file]);
			int mol_no = 0;
			for (Molecule* mol = input->read(); mol; mol = input->read(), mol_no++)
			{
				if (!mol->hasProperty(e_property))
				{
					Log.level(10) << "Compound " << mol->getName() << " in file " << names[file] << " has no score property. Skipping this compound." << endl;
					for (Size i = 0; i < new_results.size(); i++)
					{
						new_results[i]->erase(input->getCurrentLigand());
					}
					delete mol;
					continue;
				}

				double score = ((String)mol->getProperty(e_property).toString()).toFloat();

				if (score > score_cutoff || score < score_cuton)
				{
					for (Size i = 0; i < new_results.size(); i++)
					{
						new_results[i]->erase(input->getCurrentLigand());
					}
					delete mol;
//.........这里部分代码省略.........

void sortMolecules(vector<String>& names, string& output_file, Size& best_k, string& e_property, double& score_cutoff, double& score_cuton)
{
	multimap<double, Molecule*> compounds;

	for (Size file = 0; file < names.size(); file++)
	{
		GenericMolFile* input = MolFileFactory::open(names[file]);

		int mol_no = 0;
		for (Molecule* mol = input->read(); mol; mol = input->read(), mol_no++)
		{
			if (!mol->hasProperty(e_property))
			{
				Log.level(10) << "Compound " << mol->getName() << " in file " << names[file] << " has no score property. Skipping this compound." << endl;
				delete mol;
				continue;
			}

			double score = ((String)mol->getProperty(e_property).toString()).toFloat();

			if (score > score_cutoff || score < score_cuton)
			{
				delete mol;
				continue;
			}

			if ((compounds.size() < best_k || score < compounds.rbegin()->first))
			{
				compounds.insert(make_pair(score, mol));

				if (compounds.size() > best_k)
				{
					delete compounds.rbegin()->second;
					multimap<double, Molecule*>::iterator it = compounds.end();
					it--;
					compounds.erase(it);
				}
			}
			else
			{
				delete mol;
			}
			Log.level(20) << "\r" << names[file] << " : " << mol_no+1 << flush;
		}
		Log.level(20) << endl;
		Log.flush();

		input->close();
		delete input;
	}

	if (compounds.size() < best_k)
	{
		Log.level(20) << "found " << compounds.size() << " compounds matching the given criteria." << endl;
	}

	GenericMolFile* output = MolFileFactory::open(output_file, ios::out, "mol2.gz");

	for (multimap < double, Molecule* > ::iterator it = compounds.begin();
		it!=compounds.end(); it++)
	{
		*output << *it->second;
		delete it->second;
	}

	output->close();
	delete output;
}

void processMultiMoleculeFile(ScoringFunction* scoring_function, StructurePreparer* sp, String par_file, Rescoring* rescoring, bool train, double min_dock_score, String dock_score_label, GenericMolFile* input, GenericMolFile* output, double output_score_threshold, bool ignore_top)
{
	list<pair<double, bool> > rescore_list;
	double min_rescore = 1e12;
	int i = 1;

	for (Molecule* mol = input->read(); mol; mol = input->read(), i++)
	{
		String name = mol->getName();
		Log<<"====== Ligand "<<i;
		if (name != "") Log<<", "<<name;
		Log<<" ============"<<endl;
		Log.flush();

		if (mol->hasProperty("score_ligcheck"))
		{
			double score_ligcheck = ((String)mol->getProperty("score_ligcheck").toString()).toDouble();
			if (score_ligcheck < 0.95) // 0 = error, 1 = check passed
			{
				cout<<"Skipping compound because it has been marked as containing errors by LigCheck."<<endl;
				delete mol;
				continue;
			}
		}

		double score = 0;
		try
		{
			sp->prepare(mol, par_file);
			if (!rescoring || train)
			{
				scoring_function->setLigand(mol);
				scoring_function->update();
				score = scoring_function->updateScore();
				scoring_function->printResult();
			}
			if (rescoring)
			{
				if (train)
				{
					if (score > 1000)
					{
						Log.level(10)<<"score>1000, thus current compound is not added to training data set."<<endl;
						delete mol;
						continue;
					}
					rescoring->addScoreContributions(mol);
				}
				else score = rescoring->rescore(mol);
			}

			if (score < output_score_threshold)
			{
				if (ignore_top)
				{
					bool keep_dock_score = 0;
					if (mol->hasProperty("score"))
					{
						double docking_score = mol->getProperty(dock_score_label).toString().toDouble();
						if (docking_score < min_dock_score)
						{
							rescore_list.push_back(make_pair(docking_score, false));
							keep_dock_score = 1;
						}
					}
					if (!keep_dock_score)
					{
						rescore_list.push_back(make_pair(score, true));
						if (score < min_rescore)
						{
							min_rescore = score;
						}
					}
				}
				else
				{
					mol->setProperty("re-score", score);
				}
			}
		}
		catch (BALL::Exception::GeneralException e)
		{
			cout<<e.getMessage()<<endl;
			cout<<"Error! Skipping this molecule!"<<endl;
		}

		if (score < output_score_threshold && (!rescoring || !train) && !ignore_top)
		{
			*output << *mol;
		}

		delete mol;
	}

	/** If ignoring (i.e. not rescoring) the top fraction of docking-results, we need to make sure
	 that all rescored compounds get a rescore-value larger than that of the former. */
	if (ignore_top && !train)
	{
		input->reopen();
		int i = 1;
//.........这里部分代码省略.........

//.........这里部分代码省略.........
	scoring_function->update();
	scoring_function->updateScore();
	Log<<"====== Reference ligand ============"<<endl;
	scoring_function->printResult();


	/** If top fraction of docking results should not be rescored, then fetch scores and compute threshold for this fraction */

	bool ignore_top = false;
	double ignore_top_fraction = 0.0;
	if (par.has("tf"))
	{
		ignore_top_fraction = par.get("tf").toDouble();
		if (ignore_top_fraction < 1e-14 || ignore_top_fraction < 0 || ignore_top_fraction > 1)
		{
			ignore_top_fraction = 0.0;
			ignore_top = false;
		}
		else
		{
			ignore_top = true;
		}
	}
	double min_dock_score = -1e100;
	list<pair<double, bool> > rescore_list;
	if (ignore_top)
	{
		set<double> scores;
		GenericMolFile* input = MolFileFactory::open(par.get("i"));
		for (Molecule* mol = input->read(); mol; delete mol, mol = input->read())
		{
			if (mol->hasProperty("score"))
			{
				scores.insert(mol->getProperty("score").toString().toDouble());
			}
		}
		Size i = 0;
		Size max = scores.size()*ignore_top_fraction;
		set<double>::iterator s_it = scores.begin();
		for (; i < max; s_it++)
		{
			i++;
		}
		min_dock_score = *s_it;
		delete input;
	}


	/** Now, Rescore entire sd-/mol2-file   */

	double threshold = option.setDefaultReal("output_score_threshold", 1e100);

	GenericMolFile* input = MolFileFactory::open(par.get("i"));
	GenericMolFile* output = 0;

	if (simple_rescoring || !train)
	{
		output = MolFileFactory::open(par.get("o"), ios::out, input);

		DockResultFile* drf_output = dynamic_cast<DockResultFile*>(output);
		if (drf_output)
		{
			String dummy = "0";
			BALL::Docking::Result::Method method = Result::getMethod(3);
			String description = "";
			description = rescoring->getName()+"+"+scoring_function->getName();

int main(int argc, char* argv[])
{
	CommandlineParser parpars("BindingDBCleaner", "fix bindingdb.org downloads", VERSION, String(__DATE__), "Preparation");
	parpars.registerParameter("i", "input file", INFILE, true);
	parpars.registerParameter("type", "type of contained activity values: 'Ki' or 'IC50'",STRING, true);
	parpars.registerParameter("o", "output file", OUTFILE, true);
	parpars.registerParameter("target", "binding-DB target name", STRING, true);
	String manual = "This tool cleans up the sd-properties contained in sd-files downloaded from bindingdb.org.\n\nFor all compounds in the input file, the affinity value for the specified target is searched and retained but all other properties are removed. Furthermore, the IC50 or Ki value of each compound is converted to a binding-free-energy value in units of [kJ/mol] that is added as a property-tag named 'binding_free_energy'.\n\nAll compounds in the input file for which no IC50 resp. Ki value for the specified target can found, are ignored and not written to the output file.";
	parpars.setToolManual(manual);
	list<String> slist;
	slist.push_back("IC50");
	slist.push_back("Ki");
	parpars.setParameterRestrictions("type", slist);
	parpars.setSupportedFormats("i","mol2,sdf,drf");
	parpars.setSupportedFormats("o","mol2,sdf,drf");
	parpars.setOutputFormatSource("o","i");
	parpars.parse(argc, argv);

	GenericMolFile* input = MolFileFactory::open(parpars.get("i"));
	GenericMolFile* output = MolFileFactory::open(parpars.get("o"), ios::out, "mol2.gz");
	String target_name = parpars.get("target");
	target_name.trim();

	String type = parpars.get("type");
	bool use_IC50 = (type == "IC50");
	bool use_Ki = (type == "Ki");
	String response_property = "Enzymologic: IC50 nM";
	if (use_Ki) response_property = "Enzymologic: Ki nM";

	if (!use_IC50 && !use_Ki)
	{
		Log.error() << "[Error:] Please set parameter 'type' to either 'IC50' or 'Ki'!" << endl;
		return 1;
	}

	Size no_mols_found = 0;
	Size no_acitivities_found = 0;
	Size total = 1;
	Size no_invalid = 0;
	for (Molecule* mol = input->read(); mol; mol = input->read(), total++)
	{
        bool found = false;
        int prop_id = 0;
		int target_no = 1;

		for (NamedPropertyIterator it = mol->beginNamedProperty();
			it!=mol->endNamedProperty(); it++)
		{
			String prop = "TARGET Biomolecule "+String(target_no);
			if (it->getName() == prop)
			{
				String resp = response_property+" "+String(target_no);

				if (it->toString().trim() == target_name
					&& mol->hasProperty(resp))
				{
					String resp_value = String(mol->getProperty(resp).toString()).trim();

					if (resp_value == "n/a")
					{
						target_no++;
						continue;
					}

					for (NamedPropertyIterator it2 = mol->beginNamedProperty(); it2 != mol->endNamedProperty(); it2++)
					{
						it2->clear();
					}

					if (resp_value.hasPrefix(">"))
					{
						resp_value = resp_value.after(">");
					}
					if (resp_value.hasPrefix("<"))
					{
						resp_value = resp_value.after("<");
					}

					mol->setProperty("Target", target_name);
					mol->setProperty(response_property, resp_value);
					double value = resp_value.toFloat();
					if (use_IC50)
					{
						String name = "Enzymologic: pIC50 nM";
                        if (prop_id > 0)
                        {
                            name += " "+String(prop_id);
                        }
						mol->setProperty(name, -log10(value));
					}
					else if (use_Ki)
					{
						String name = "Enzymologic: pKi nM";
                        if (prop_id > 0)
                        {
                            name += " " + String(prop_id);
                        }
						mol->setProperty(name, -log10(value));
					}

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

//.........这里部分代码省略.........
    {
      // Option 1: Generate output file names from specified pattern

      String pattern = parpars.get("outname_pattern");
      
      for (unsigned int i=0; i!= n_outfiles; ++i)
      {
        outfile_names.push_back(getOutputFileName(pattern, true, infile_type, i));
      }
    }
    else
    {
      // Option 2: Simple indexing of input file name
      
      for (unsigned int i=0; i!= n_outfiles; ++i)
      {
        outfile_names.push_back(getOutputFileName(infile_name, false, infile_type, i));
      }
    }
      
  }
  
  
  // Now do the splitting
  
  input = MolFileFactory::open(infile);

  drf_input = NULL;
	drf_input = dynamic_cast<DockResultFile*>(input);
	if (drf_input) 
  {
    drf_input->selectAllResultsForInput();
  }
  
  for (unsigned int i=0; i!=outfile_names.size(); ++i)
  {
    conformation_ids.clear();
    
    output = MolFileFactory::open(outfile_names[i], ios::out, infile_type);
    drf_output = dynamic_cast<DockResultFile*>(output);
    if (drf_input && drf_output)
    {
      drf_output->disableAutomaticResultCreation();
    }

    
    if (i < outfile_names.size()-1)
    {
      // Not the last file - so number of molecules is standard
      
      for (unsigned int j=0; j!=mpf; ++j)
      {
        mol = input->read();
        if (drf_input && drf_output && mol->hasProperty("Conformation_input_UID"))
        {
          conformation_ids.insert(mol->getProperty("Conformation_input_UID").toString());
        }

        output->write(*mol);
        delete mol;
      }
    }
    else
    {
      // Last output file - so write remaining molecules into it
      
      mol = input->read();
      while (mol)
      {
        if (drf_input && drf_output && mol->hasProperty("Conformation_input_UID"))
        {
          conformation_ids.insert(mol->getProperty("Conformation_input_UID").toString());
        }

        output->write(*mol);
        delete mol;
        mol = input->read();
      }
    }

    if (drf_input && drf_output)
    {
      const vector<Result*>* results = drf_input->getResults();
      for (unsigned int i = 0; i < results->size(); ++i)
      {
        Result* new_res = new Result(*(*results)[i], conformation_ids);
        drf_output->addResult(new_res);
      }
    }
    
    output->close();
    delete output;
  }  

  input->close();
  delete input;
  
  
  return 0;
}