From 4988549536659be427261a9bdbf56084a4352fe7 Mon Sep 17 00:00:00 2001
From: Baptiste LENIAU <baptiste.leniau@subatech.in2p3.fr>
Date: Thu, 26 Jun 2014 13:13:36 +0000
Subject: [PATCH] Add GetMass and MeanMolarMass Function

git-svn-id: svn+ssh://svn.in2p3.fr/class@276 0e7d625b-0364-4367-a6be-d5be4a48d228
---
 .../branches/CLASSV3/src/IsotopicVector.cxx   | 45 ++++++++++++++++++-
 source/branches/CLASSV3/src/ZAI.cxx           | 22 ++++-----
 source/branches/CLASSV3/src/ZAIMass.cxx       | 14 +++++-
 3 files changed, 69 insertions(+), 12 deletions(-)

diff --git a/source/branches/CLASSV3/src/IsotopicVector.cxx b/source/branches/CLASSV3/src/IsotopicVector.cxx
index 0bb82cc29..8441ad37d 100755
--- a/source/branches/CLASSV3/src/IsotopicVector.cxx
+++ b/source/branches/CLASSV3/src/IsotopicVector.cxx
@@ -8,7 +8,7 @@
 #include <fstream>
 #include <string>
 #include <algorithm>
-
+#include "CLASSHeaders.hxx"
 	//________________________________________________________________________
 	//________________________________________________________________________
 	//
@@ -451,6 +451,49 @@ IsotopicVector	IsotopicVector::GetSpeciesComposition(int z) const
 	return IV;
 	
 }
+//________________________________________________________________________
+double IsotopicVector::GetTotalMass() const
+{
+	double TotalMass = 0;
+	double AVOGADRO=6.02214129e23;
+	map<ZAI ,double >::iterator it;
+	map<ZAI ,double > isotopicquantity = GetIsotopicQuantity();
+	for( it = isotopicquantity.begin(); it != isotopicquantity.end(); it++)
+	{
+		ZAI zai = ((*it).first);
+		double MolarMass = zai.GetMass();
+		TotalMass += (*it).second/AVOGADRO * MolarMass ;
+	}	
+		
+
+
+	return TotalMass*1e-6;//in tons
+
+}
+//________________________________________________________________________
+
+double IsotopicVector::MeanMolar() const
+{
+	double MeanMolar = 0;
+	map<ZAI ,double >::iterator it;
+	map<ZAI ,double > isotopicquantity = GetIsotopicQuantity();
+
+	double NTot=0;
+	for( it = isotopicquantity.begin(); it != isotopicquantity.end(); it++)
+		NTot+= (*it).second ;
+
+	for( it = isotopicquantity.begin(); it != isotopicquantity.end(); it++)
+	{
+		ZAI zai = ((*it).first);
+		double MolarMass = zai.GetMass();
+		MeanMolar += (*it).second/NTot * MolarMass ;
+	}	
+		
+
+	return MeanMolar;
+
+}
+//________________________________________________________________________
 
 vector<ZAI> IsotopicVector::GetZAIList() const
 {
diff --git a/source/branches/CLASSV3/src/ZAI.cxx b/source/branches/CLASSV3/src/ZAI.cxx
index 4f795aa55..86e81c5b9 100755
--- a/source/branches/CLASSV3/src/ZAI.cxx
+++ b/source/branches/CLASSV3/src/ZAI.cxx
@@ -1,4 +1,6 @@
 #include "ZAI.hxx"
+#include "CLASSHeaders.hxx"
+
 
 //const string DEFAULTDATABASE = "DecayBase.dat";
 //________________________________________________________________________
@@ -18,6 +20,7 @@ ZAI ZAI::operator=(ZAI IVa)
 	fZ = IVa.Z();
 	fA = IVa.A();
 	fI = IVa.I();
+	fMass =  IVa.GetMass();
 	return *this;
 }
 
@@ -39,21 +42,20 @@ ZAI::~ZAI()
 }
 
 //________________________________________________________________________
-ZAI::ZAI(int Z, int A, int I)
+ZAI::ZAI(int Z, int A, int I,bool IsMassSet)
 {
 		
 	fZ=Z;
 	fA=A;
 	fI=I;
-		
-}
 
-//________________________________________________________________________
-double ZAI::GetMass()
-{
-		
-	return fMass;
-}
+	if(!IsMassSet)
+		fMass=0;
 
-//________________________________________________________________________
+	else
+	{
+		fMass=cZAIMass.GetMass(Z,A);
+	}
+
+}
 
diff --git a/source/branches/CLASSV3/src/ZAIMass.cxx b/source/branches/CLASSV3/src/ZAIMass.cxx
index e44d880af..29be907af 100644
--- a/source/branches/CLASSV3/src/ZAIMass.cxx
+++ b/source/branches/CLASSV3/src/ZAIMass.cxx
@@ -29,7 +29,7 @@ ZAIMass::ZAIMass()
 	while (infile>>Z>>A>>Name>>MassUnity>>MassDec>>error)
 	{
 		double Masse=MassUnity+MassDec*1e-6;
-		fZAIMass.insert( pair< ZAI,double >( ZAI(Z,A,0), Masse ) );
+		fZAIMass.insert( pair< ZAI,double >( ZAI(Z,A,0,false), Masse ) );
 	}
 
 	infile.close();
@@ -42,3 +42,15 @@ ZAIMass::~ZAIMass()
 	fZAIMass.clear();
 }
 
+
+double ZAIMass::GetMass(const int Z,const int A) const
+{
+	map<ZAI,double> ZAIMasscpy = fZAIMass ;
+	map<ZAI,double>::iterator  MassIT = ZAIMasscpy.find( ZAI(Z,A,0,false) );
+	if(MassIT==fZAIMass.end())
+		return A;
+	
+	else
+	   return MassIT->second;
+
+}
-- 
GitLab