The code for G4NeutronHPContAngularPar::Sample() has a bug located in the for loop at line 550 in G4NeutronHPContAngularPar.cc (version 4.9.5p01).  

This code is only used when dealing with outgoing angles that are saved in the G4NDL format as a table of cos(theta) values, with corresponding probabilities, that are indexed by the incoming neutron energy.  Using the default G4NDL data (G4NDL 3.14 or 4.0), this is not a problem because the outgoing angle is described by legendre polynomials.  However, when using different neutron data, the following code causes a runtime error:


	for(i=0; i<nAngularParameters; i++)
	{
		theBuff1.SetX(i, theAngular[it-1].GetValue(i));
		theBuff1.SetY(i, theAngular[it-1].GetValue(i+1));
		theBuff2.SetX(i, theAngular[it].GetValue(i));
		theBuff2.SetY(i, theAngular[it].GetValue(i+1));
		i++;
	}


The first problem is that "i" increments by two for every iteration of the loop.  Therefore, the SetX() and SetY() functions try to set elements [0,2,4,...] of the buffer objects "theBuff1" and "theBuff2".  This causes the runtime error since elements are skipped in these buffers.  Additionally, the angular data, "theAngular", has the following format

	theAngular = {P_{energy},cos_1(theta),P_1,cos_2(theta),P_2,...}

where "P_{energy}" is the probability of the outgoing neutron energy, "cos_i(theta)" is the ith outgoing angle, and "P_i" is the probability of the ith outgoing angle.  Therefore, the GetValue() function should start at 1, not 0.

Fixing this bug is relatively simple, and the solution I used is shown below


	for(i=0, j=1; i<nAngularParameters; i++, j+=2)
	{
		theBuff1.SetX(i, theAngular[it-1].GetValue(j));
		theBuff1.SetY(i, theAngular[it-1].GetValue(j+1));
		theBuff2.SetX(i, theAngular[it].GetValue(j));
		theBuff2.SetY(i, theAngular[it].GetValue(j+1));
	}


In the fixed code, a separate index "j" is used to get the cos(theta) data from theAngular.