MODULE input
	IMPLICIT NONE
	INTEGER::samples,cell_XY,cell_Z,min_cell_E,max_cell_E,min_A,max_A,min_Z,max_Z
	REAL(8)::r_galaxy,z_galaxy,grid,energy_grid,D0,delta,per_disc_type1a,per_disc_type2,energy_start,energy_end
	REAL(8)::inner_r_4arm_set,pitch_angle_4arm_set,width_4arm_set,omega_4arm_set,start_angle_4arm_set,per_4arm_set
	REAL(8)::inner_r_2arm_set,pitch_angle_2arm_set,width_2arm_set,omega_2arm_set,start_angle_2arm_set,per_2arm_set
	REAL(8)::per_Ni,Per_Fe,per_Ti44,per_Fe60,per_Si,per_Mg,per_Ne,per_O,per_N,per_C,per_e,per_P
	REAL(8)::half_t_Ni59,half_t_Fe60,half_t_Mn53,half_t_Al26,half_t_Be10,r_plot,ECdelta,ECa
	CHARACTER::model*11,calculate*6,plot*9,density_plot*10,crossection*6
	REAL(8), PARAMETER::megayear_sec=3.16*1d13,v_light_v_cm=3*1d10,cm2_mb=1d27,erg_GeV=624.150974
	REAL(8), PARAMETER::pi=3.1415926536,r_e=2.818*10**(-13.),m_e=9.1*10**(-28.),kpc_km=30.857*1d15,kpc_cm=30.857*1d20
	
	CONTAINS
	SUBROUTINE setup()
		IMPLICIT NONE
		REAL(8)::age_4arm_set,density_energy
		samples=1d8
		model='spiral arms'						!disc, spiral arms
		calculate='S/P'							!S/P, galaxy
		plot='PAMELA'								!B/C, Iron, PAMELA, nuclei, electrons, protons
		density_plot='anisotropy'					!X-Y, X-Z, Z-phi, gamma, anisotropy
		crossection='YIELDX'					!YIELDX, WEBBER
		SELECT CASE(model)
			CASE('disc')
				SELECT CASE(plot)
					CASE('B/C')
						delta=0.4				!nominal value:0.4 disc B/C
						D0=8*1d27				!nominal value:8*1d27 disc B/C
					CASE('Iron')
						SELECT CASE(crossection)
							CASE('YIELDX')
								delta=0.5		!nominal value:0.5 YIELDX disc Iron
								D0=4*1d27		!nominal value:4*1d27 YIELDX disc Iron
								ECa=0.000038
								ECdelta=1.8
							CASE('WEBBER')
								delta=0.4		!nominal value:0.4 WEBBER disc Iron
								D0=7*1d27		!nominal value:7*1d27 WEBBER disc Iron
								ECa=0.000038
								ECdelta=1.8
						END SELECT
					CASE DEFAULT
						delta=0.45
						D0=6*1d27
				END SELECT
				z_galaxy=1						!nominal value:1 disc
			CASE('spiral arms')
				delta=0.4						!nominal value:0.4 spiral arms
				D0=1.2*1d27						!nominal value:1.2*1d27 spiral arms
				z_galaxy=0.25					!nominal value:0.25 spiral arms
				age_4arm_set=5					!nominal value:age_4arm_set=5
				ECa=0.000068*((z_galaxy/0.5)**(0.24))*((age_4arm_set/10)**(-0.28))	!nominal value:normalization=0.000068,z_power=0.24,age_power=-0.28
				ECdelta=1.8						!nominal value:ECdelta=1.8
		END SELECT
		r_galaxy=15								!nominal value:r_galaxy=15
		D0=(D0/kpc_cm**2)*megayear_sec
		SELECT CASE(model)
			CASE('disc')
				per_disc_type1a=19.3
				per_disc_type2=80.7
				per_disc_type2=per_disc_type2/100
				per_disc_type1a=per_disc_type1a/100
			CASE('spiral arms')
				pitch_angle_4arm_set=28
				pitch_angle_2arm_set=11
				per_disc_type1a=0!19.3
				per_disc_type2=0!8.1
				per_2arm_set=0!24.2
				per_4arm_set=100!48.4
				width_2arm_set=0.1!6
				width_4arm_set=0.1!8
				omega_2arm_set=25
				omega_4arm_set=15
				inner_r_2arm_set=4
				inner_r_4arm_set=6.5
				per_disc_type1a=per_disc_type1a/100
				per_disc_type2=per_disc_type2/100
				per_2arm_set=per_2arm_set/100
				per_4arm_set=per_4arm_set/100
				pitch_angle_2arm_set=pitch_angle_2arm_set*pi/180
				pitch_angle_4arm_set=pitch_angle_4arm_set*pi/180
				start_angle_2arm_set=-3.234
				start_angle_4arm_set=-(TAN(pi/2-pitch_angle_4arm_set)*LOG(8.5/inner_r_4arm_set)-age_4arm_set*3*pi/(5*180))
		END SELECT 
		SELECT CASE(plot)
			CASE('B/C')
				per_Si=2				!0% disc, 2% spiral arms
				per_Mg=2				!0% disc, 2% spiral arms
				per_Ne=4				!0% disc, 4% spiral arms
				per_O=51				!55% disc, 51% spiral arms
				per_N=4					!5% disc, 2% spiral arms
				per_C=37				!40% disc, 37% spiral arms
				per_Si=per_Si/100
				per_Mg=per_Mg/100
				per_Ne=per_Ne/100
				per_O=per_O/100
				per_N=per_N/100
				per_C=per_C/100
				half_t_Al26=0.7166
				half_t_Be10=1.387
				half_t_Al26=half_t_Be10/LOG(2.)
				half_t_Be10=half_t_Be10/LOG(2.)
				min_A=9
				max_A=30
				min_Z=4
				max_Z=14
			CASE('Iron')
				per_Ni=5				!5% disc, 5% spiral arms
				per_Fe60=0.0035				!94.8% disc, 94.8% spiral arms
				per_Ti44=0.38			!0.38% disc, 0.38% spiral arms
				per_Fe=95-(per_Ti44+per_Fe60)
				per_Ni=per_Ni/100
				per_Fe=per_Fe/100
				per_Fe60=per_Fe60/100
				per_Ti44=per_Ti44/100
				half_t_Ni59=0.076
				half_t_Fe60=2.6
				half_t_Mn53=3.7
				half_t_Ni59=half_t_Ni59/LOG(2.)
				half_t_Fe60=half_t_Fe60/LOG(2.)
				half_t_Mn53=half_t_Mn53/LOG(2.)
				min_A=44
				max_A=64
				min_Z=20
				max_Z=29
			CASE('PAMELA')
				per_e=0.7
				per_P=99.3
				per_e=per_e/100
				per_P=per_P/100
				min_A=0
				max_A=1
				min_Z=-1
				max_Z=1
		END SELECT 
		SELECT CASE(calculate)
			CASE('S/P')
				energy_start=1
				energy_end=100000
			CASE('galaxy')
				grid=0.1
				cell_XY=(r_galaxy/grid)+1
				cell_Z=(z_galaxy/(grid/10))+1
				SELECT CASE(density_plot)
					CASE('anisotropy')
						energy_start=1
						energy_end=100000
					CASE default
						density_energy=1
						energy_start=density_energy

				END SELECT 
		END SELECT 
		energy_grid=0.1
		IF(energy_start<1)THEN
			min_cell_E=(LOG10(energy_start)/energy_grid)-1
		ELSE IF(energy_start>1)THEN
			min_cell_E=(LOG10(energy_start)/energy_grid)+1
		ELSE
			min_cell_E=(LOG10(energy_start)/energy_grid)
		END IF
		SELECT CASE(calculate)
			CASE('S/P')
				IF(energy_end<1)THEN
					max_cell_E=(LOG10(energy_end)/energy_grid)-1
				ELSE IF(energy_end>1)THEN
					max_cell_E=(LOG10(energy_end)/energy_grid)+1
				ELSE
					max_cell_E=(LOG10(energy_end)/energy_grid)
				END IF
			CASE('galaxy')
				SELECT CASE(density_plot)
					CASE('anisotropy')
						IF(energy_end<1)THEN
							max_cell_E=(LOG10(energy_end)/energy_grid)-1
						ELSE IF(energy_end>1)THEN
							max_cell_E=(LOG10(energy_end)/energy_grid)+1
						ELSE
							max_cell_E=(LOG10(energy_end)/energy_grid)
						END IF
					CASE default
						max_cell_E=min_cell_E+3
						r_plot=8.5
				END SELECT 
		END SELECT 
	END SUBROUTINE setup
END MODULE input
