/* 
  lunar-ice-cor2.c, contains comprehensive designs of lunar rovers...

   Possible Changes:
   -> temp variables in functions?
   -> redundant distrib/distribution syntax !
   -> UI ability to change label - element, component, unit, etc
   -> Function overloading?
   -> When it is 'key mass', the origin should default to the current component,
   not root
   -> add keyword parent, as a reference, eg key mass of parent
   -> ability to key off grep expressions eg. sum(key boog* /suboog* of ...)
   so that people can operate on name groups rather than on just arrays
   -> Is type automatically inherited? - NO
   -> can a default be defined as a formula? say a function of other defaults?
   -> automate mass, power, ... calc given N, no of such components
   (define keyword value and allow defaults to be defined in terms
   of user defined "ivalue"), define function ivalue() to access that
   ivalue as well. so mass would become mass = real ivalue*N:0 kg;
   ivalue functions apply only when it is defined, reverts to default
   otherwise. dvalue would refer to derived value, and would be
   the one accessed by default.
   -> Probably need default consistency constraints as well.
   Maybe a default component ought to be defined which is included
   in all components. Ability to "fill in" values that are undefined
   from others which are, using default relations is useful.
   For e.g. given density, mass and two dimensions, third would be obvious
   -> To do the above an ideal thing is an equality constraint, and a symbolic
   algenra operator which automatically derives the missing values...
   -> need a zero finding function - zero(fn,init,a,b) [a,b] optional range
   -> Forgive plural "s" in units. Consequently disallow xyz & xyzs definitions
   -> Does adding lots of defaults add inefficiency to sums (0 defaults)
   and products (1 default). Also confirm that a global sum is not effected
   w/wo default
   -> Optimization, if a calculated quantity is not used by others. Do not
   calc until asked for.
   -> Need a "lookup" function which takes two arrays of numbers, an interpolation
   function tag (nearest, linear, quad, spline), a variable in one of the axes; and
   returns the corresponding value from the other axis. Could be used for battery
   specs.
   -> Once compiled, fold constants, for eg. if n = 3.2, n becomes a constant
   wherever used, tag it as such, any expression is a constant until a variable
   is used in it, and so on...

   -> Ability to incrementally ReCompile is essential for speed && felxibility

   Is it possible to switch pointers/ templates using if()?
   Is it OK to use formula in initial value? NO
   Does modify or template inclusion allow type def?
   
   Are C, K defined as centigrade and Kelvin?
   Add automated computing budgets, software modules and computing needs.
   Need automated comm budget
   What is the Syntax for type casting?
   Does transparent make it invisible to access or just set operations?
   
   type syntax
   type inherited A, local B; // local is default
   
   Get Mf/Mi from Gordy, RTG info from AMPS

   Freeze All Options to improve speed

   The pointer can be reassigned once. Explain, for e.g. compression as a
     function of frequency?

     Problem - Latency needs typecasting, does not fall out as a formula
       - find an elegant representation if u can

  4436, 4380
   */

/**** OPTIMIZATION POSSIBILITIES
  Locomotion - d, w, n
  Comm, P, Dia/N; parabolic vs phased array (play w tiltmax)
  Battery type
  
  
************/

default mass = real 0 kg;
default power = real 0 watt;
default volume = real 0 m3;
default T = real 0 kelvin;
default cost = 0 megadollar;
default dcost = 0 megadollar;
default efficiency = real 1;
default duration = real 0 sec;
default speed = real 0 [m/s];
default freq = real 0 hz;
default theta = real 0 radian;
default N = int 1;
/* default mips* = 0 byte/s; // is this valid ? */
default mipsGint = 0 ips;
default mipsGfp = 0 ips;
default mipsDint = 0 ips;
default mipsDfp = 0 ips; /* does this cause inefficiency? */
default resolution = 0 byte;
default downlink = 0 bps;
default reliability = 0.99; /* over mission lifetime */
default latency = real 0.0 sec;
/* no of components, need to automate mass, power, volume, etc ? */
default area = 0 m2;

Pi = 3.14159265358979323846 radian;
triple y_axis = triple( 0, 1, 0 );

Component LunarRover {
  
  Component Mission : Icediscovery {} 
  Component Launcher : DeltaII7925 {} 
  
  component Earth
    {
      /* this would normally be added by the robot design rather
	 than in the component itself. Do when we have ability
         to reference & include components */
      
      g = 9.81 [m/s2];
      
      CommGainAtmosphere  = -0.3;
      CommGainIonosphere  = -0.1;
      distrib CommGainRain = distribution(max=0, min=-1, avg=-0.2);
    }
  
  Component Light
    {
      speed = 2.998E8 [m/s];
    }
  
  Component Planet : Moon{}
  
  component Space
    {
      T = 3 kelvin; /* degree Kelvin */
    }
  
  Component LauncherInterface {mass = 0.005*(old spacecraft.mass);
			       cost = 0.2;}
      
  Component Spacecraft
    {
      transparent mass = sum( key mass of Spacecraft ) 
	+ (Lander.LanderRoverMass - Lander.mass - Rover.mass):
      1500 kg, clock coriolis, main_clock;
      H = SRM.H + Lander.H + Rover.H;
      CGz = (SRM.mass * SRM.CGz + Lander.LanderRovermass * 
	     (SRM.H + Lander.LanderRoverCGz))/
	       (SRM.Mass + Lander.LanderRoverMass);
      
      Launcher.massTLI <= Spacecraft.mass + LauncherInterface.mass;
      
      Component Rover
	{
	  transparent mass = sum(key mass of Rover) : 392 kg;
	  /* Bug below ?, dependency not auto found */
	  weight = old mass * planet.g;
	  transparent power = sum(key power of Rover) : 400 watt; /* W */
	  transparent Nactuators = sum(key Nactuators of Rover);
	  transparent Nsensors  = sum(key N of type sensor of Rover);
	  L = Locomotion.L;
	  W = Locomotion.W;
	  H = wheel.d+ebox.h+mast.h;
	  CGmass = locomotion.mass + ebox.mass + mast.mass + mast.payload.mass;
	  
	  CGz = (locomotion.mass*wheel.d/2 + ebox.mass*(wheel.d+ebox.h/2)
		 + mast.mass*(wheel.d+ebox.h+mast.h/2) + mast.payload.mass*
		 (wheel.d+ebox.h+mast.h) +
		 (old rover.mass-CGmass)*(wheel.d+ebox.h))/
		   (old rover.mass); // from rover bottom
	  /* meters offset from bottom center */
	  CGLoff = 0.1*Locomotion.L; 
	  CGwoff = 0.1*Locomotion.W;
	  /* length & width wise horizontal offsets */

	  Component Locomotion
	    {
	      Component Soil  : LunarSoil {}
	      
	      
	      /* these, wheel.d,w have to change for small rovers ? */
	      Component TempWheel
		{
		  no = 4; /* n wheels */
		  OFactor = if(no==4, 0.25, 0.30);
		  dia = max(mission.obstaclemax/OFactor, wheel.sinkMax*4.0, 
			    old wid) : 0.5 m;
		  wid = max(wheel.pressuremax/soil.pressuremax * old wid, dia*0.1) : 0.2 m;
		  /* May be worth checking out flexible wheel for power/ size */
		  Component Wheel  : RigidWheel {}
		}
	      Lnew = max( dia*n / 2 + 
			  dia*( -1 + n / 2 )*0.3, 
			 L * slopeLmax/mission.slopemax )
         	: 2.0 m;
	      L = old Lnew; /* length */
	      Wnew = max(wid*3, W * slopeWmax/mission.slopemax) : 2.0 m; /* width */
	      W = old Wnew; /* width */
	      Nactuators = wheel.n + 2; // ?
	      
	      slopeLmax = atan((L/2 - wheel.d/2 - rover.CGLoff)/rover.CGz);
	      slopeWmax = atan((W/2 - wheel.w*0.1 - rover.CGwoff)/rover.CGz);
	      slopemax = max(slopeLmax, slopeWmax);
	      
	      weight = (old Rover.mass )*planet.G;
	      wheelLoadMax = weight/2;
	      
	      /* Total electric Driving power in W */
	      
	      P(theta=0 radian,V=mission.speed) = amplifier.P(theta,V)*wheel.n;
	      
	      /* Now calc electric power in for different conditions ? */
	      /* Drive sizing to include flat and slope conditions ? */
	      
	      mass = 0.30/(1-0.30)*(old Rover.mass);
	      power = P(theta=0,V=get_avg(mission.speed))*1.05 : 90 Watt; 
	      /* 1.05 is margin */
	      cost = 6;
	      
	      Component Axle
		{
		  type wheel;
		  efficiency = 0.8;
		}
	      Component Bearing
		{
		  type wheel;
		  efficiency = 0.9;
		}
	      Component SteerLink{}
	      
	      
	      Component Amplifier
		{
		  type inEBox;
		  /* electrical power input to amplifier */ /* Nm */
		  N = locomotion.Nactuators;
		  P(theta,V) = Wheel.Tw(slope = theta)*
		    wheeldrive.omega(V)/efficiency/wheeldrive.Nu;
		  
		  /* At peak torque condition */
		  PTmax = wheel.Torquemax*wheeldrive.omegapeak/efficiency/Wheeldrive.Nu;
		  Volume = PTmax / ( 1000 [Watt/liter] ) * N;
		  efficiency = 0.85;
		}
	      Component WheelDrive
		{
		  type inherited wheel;
		  
		  Component reduction1
		    {
		      efficiency = 0.7; /* 0.7 */
		      ReductionRatio = 100; /* Reduction ratio */
		    }
		  Component reduction2
		    {
		      efficiency = 0.7; /* 0.7 */
		      ReductionRatio = 100; /* Reduction ratio */
		    }
		  safetyFactor = 1.5; /* for specing torque */
		  
		  /* Overall efficiency of wheel & drive train in delivering
		     power from motor output to terrain */
		  Nu = product(key efficiency of type wheel of locomotion);
		  
		  /* Nominal drive train spec at output */
		  Torque = wheel.Tw();
		  omega(V=mission.speed) = V/(wheel.D/2-wheel.defl);
		  
		  /* Peak drive train spec at output */
		  TorquePeak = wheel.Torquemax;
		  omegaPeak = get_min(mission.speed)/(wheel.D/2-wheel.defl);
		  
		  /* software */
		  safeOpsFP = 200 ips * omega(mission.speed)/(1 Hz);        

		  /* CHANGED to typecast (ips) omega does not work ?  
		   omega does not work, changed to omega(), please fix it ?*/
		  
		  controlINT = 100 * omega()/(1 Hz) * 25 ips; /* every 2 deg, control loop */
		  /* CHANGED next 2 lines */
		  MipsGint = ControlInt * Nactuators; /* actuator control loop */
		  MipsGfp = safeOpsFP * Nactuators;
		  
		  /* constraints */
		  MotorNcycles = 10E6; /* real numbers needed */
		  MotorNCycles > wheel.Ncycles*2.0*ReductionRatio;
		} /* WheelDrive */
	      
	      Component SteerDrive
		{
		  efficiency=0.3;
		}
	      
	      /* software */
	      safeingFP = 500 ips* mission.speed/(0.05 [m/s]); /* check every 2 cm */
	      controlFP = 1000 ips* mission.speed/(0.01 [m/s]);
	      /* coordinated control every 2 cm */
	      MipsGfp = safeingFP + controlFP; /* CHANGED to typecast */
	      
	      /* constraints */
	      sqrt(L^2 + W^2) <= launcher.Dia;
	      wheel.N/2.0*1.1*wheel.D <= L;
	      slopemax >= mission.slopemax;
	      /* the step/D ratio is really a function of (L+-CGLoff)/L,
		 traction/weight & L/D, get something better for the hack */
	      D >= mission.obstacleMax/0.25; /* What is this factor really? */
	    } /* locomotion */
	  
	  Component Structure
	    {
	      mass = 0.2*(old Rover.mass
			  /* except desc of locomotion ? */
			  - locomotion.mass) /* HACK ?.  need to except self */ ;
	      cost = 2;
	      Component Body{}
	      Component EBox
		{
		  type Insulated;
		  
		  transparent mass = sum(key mass of type inEBox of Rover);
		  /* volume solid is volume of contents */
		  volumesolid = sum(key volume of type inEBox of Rover);
		  transparent volume = volumesolid * (1.0/0.4 - 1.0); /* CHANGED */
		  /* 0.4 is packing fraction */
		  H = ((volume + volumesolid)/(1 m3))^(1./3.0)*1.0 m : 1.0 m; 
		  /* BUG? will the above work for units, hacked */
		  Area = 6.0*(old H)^2; // assuming it is a cube
		  
		  MipsGfp = 100 * 1 ips; /* one safety check & log per second */
		  
		  /* Max temp EBox can take is min of all Tmax'es */
		  Tmax = Min(key Tmax of type inEBox of Rover);
		  /* similarly for Tmin */
		  Tmin = Max(key Tmin of type inEbox of Rover);
		}
	    } /* Structure */
	  
	  Component PowerSystem
	    {
	      transparent mass = sum(key mass of Powersystem);
	      powerNeed = old Rover.power;
	      Component PowerSource{} //  : power.SolarArray {}
	      
	      Component EnergyStorage
		{
		  type inherited inEBox;
		  EnergyNeed = mission.energyNeed;
		  
		  /* add equation for battery life using log scale ? */
		  Component BatteryPack  : AgZnHR190 {}
		}
	      
	      Component PMAD
		{
		  type inherited inEBox;  
		  /* Power Management And Distribution */
		  cost = 1;
		  margin = 0.2; /* power margin for the overall system */
		  /* Mass of a regulator of given power, assume factor of safety 2 */
		  
		  Component Regulator {
		    mass = (2* old rover.power/1000 watt) * 4 kg;
		    volume = (mass/4 kg)*(30 cm * 30 cm * 30 cm);
		  }
		  Component DCDC
		    {
		      efficiency = 0.85;
		      mass = dcpowerneed/(15 watt)* 50 gram;
		      volume = (mass/50 gram)*(3 cm * 3 cm * 1 cm);
		      dcpowerneed = old Rover.power - old locomotion.power
			- old power;
		      /* actually, need to seperate dc-dc power from margin,
			 num correct though ? */
		      power = dcpowerneed/efficiency * (1+margin) -
			dcpowerneed + locomotion.power*margin : 50 watt;
		      /* this is PMAD power */
		    }
		  Component PowerBus
		    {
		      mass = (old EBox.H/1.5 m)*(dcPowerNeed/300 watt)*3 kg;
		      volume = old EBox.H * 10 cm * 10 cm;
		    }
		  
		  MipsGfp = 1000 ips;
		  /* one safety check & adjustment per second,
		     any spikes, transients are handled by regulator/ safety circuit.
		     any swtich operation is a triggered operation */
		  
		}/* PMAD */
	    } /* PowerSystem */
	  
	  Component Thermal
	    {
	      transparent mass = sum(key mass of thermal);
	      cost = 2;
	      Component MLIblanket
		{
		  transparent area = sum(key area of type Insulated of Rover);
		  Nlayers = 10;
		  Component insulation  : thermal.MLI {}
		}
	      
	      Component TempHeatPipe
		{
		  type inEBox, insulated;
		  H = old EBox.H;
		  Component HeatPipe  : ConvectiveHeatPipe {}
		}
	      Component EBoxRadiator
		{
		  dayHeat = sum(key power of type inEBox of Rover) -
		    insulation.heatflow(Tebox,mission.Tmax, Ebox.area);
		  dayH1 = sum(key power of type inEBox of Rover);
		  dayH2 = insulation.heatflow(Tebox,mission.Tmax, Ebox.area);
		  /* inflow of heat */
		  Tebox = EBox.Tmax;
		  Thot  =  Tebox + heatpipe.dT(dayHeat); // radiator hot temp
		  Tcold = EBox.Tmin;
		  
		  
		  Component radplate  : Radiator {} /* Radiator plate */
		  
		  Area = radplate.sarea(Thot,space.T,dayHeat, mission.latitude);
		  heatflow(Thot,space.T, area) >= dayHeat;
		  
		}
	      
	      nightPowerNeed = if(EBox.Tmin < mission.Tmin,
				  (radplate.heatcold(Tcold,space.T, EBoxRadiator.area) +
				   insulation.heatflow(Tcold,space.T))*2.0,
				  0.0 watt);
	      /* 2.0 is factor of safety */
	      
	      nightPowerNeed <= sum(key NightPowerOut of rover);
	      Component HeatSource  : HeatSourceRHU {} // Or HeatSourceNone
	      /* if this is not enough, battery takes the slack */
	      
	      mipsGfp = 200 * 1 ips; /* one safety check and log per second */
	    } /* Thermal */
	  
	  Component Computronics
	    {
	      type inherited inEBox;

	      transparent mass = sum(key mass of Computronics);
	      
	      Component Cage{
		mass = boards.nboards/6*4 kg;
		volume = boards.nboards*11 inch * 9 inch * 1 inch; // Correct 6 U?
		cost = 0.3;}
	      
	      /* once there ability to spawn or prototype boards,
		 spawn N boards, and tweak power on each to on/off ? */
	      Component Boards{
		
		OSgOverhead = 2.0;
		OSdOverhead = 1.3;
		mipsGIntTotal = sum(key mipsGInt of rover) * (OSgOverhead-1);
		/* General purpose needed */
		mipsDIntTotal = sum(key mipsDInt of rover) * (OSdOverhead-1);
		/* Other mips, DSP or otherwise */
		mipsGfpTotal  = sum(key mipsGfp of rover) * OSgOverhead;
		mipsDfpTotal  = sum(key mipsDfp of rover) * OSdOverhead;
		
		component messaging { /* inter,intraboard messaging */                                                 
		  freq = 100 Hz; 
		  /* each board is hit with 10 mips overhead, OK?  */
		  mipsGint = freq * old boardsactive * 100000 byte;}
		
		/* use ceil, instead of + 1? */
		int boardsActive = 
		  min((mipsGIntTotal/bmipsGint + mipsGFpTotal/bmipsGfp)
		      /efficiencyG+1,
		      (mipsDIntTotal+mipsGIntTotal)/
		      (bmipsDint*efficiencyD + bmipsGint*efficiencyG)
		      + (mipsDFpTotal+mipsGFpTotal)/
		      (bmipsDFp*efficiencyD + bmipsGFp*efficiencyG) + 1) : 1;
		
		/* ?Hack, ideally same should be done for Nactuators (D/A),
		   however, usually shortage is only for AtoD */
		SensorsNeeded = rover.nsensors - boardsactive*board.NAtoD;
		int Nsensorboards = 
		  if(SensorsNeeded > 0, SensorsNeeded/8 + 1,0);
		sbPower = 8 Watt;
		
		int Nboards = boardsActive+1 + Nsensorboards;
		mass = bmass*N; /* including heat sink, 5 of these */
		power = bpower*boardsActive + NsensorBoards*SBpower;
		
		/* lots of $ to develop, mainly for radiation readiness */
		cost = board.costNRE + N*0.15 megadollar;
		
		
		Component Board 
		  {
		    Component CPU  : mipsR10000 {}
		    Component DSP  : TIC60 {}
		    /* old board use R3000, empty for DSP */
		    NAtoD = 16;
		    NDtoA = 4;
		    bmass = 1 kg; /* including heat sink */
		    bpower = CPU.bpower + DSP.bpower + 5 watt + 5 watt;
		    /* 5W for memory, 5 W for I/O */
		    costNRE = CPU.dcost + DSP.dcost + 0.5 megadollar;
		    /* in addition to chip cost, board dev +
		       mil-specization - temp, vibration + vacuum */
		    
		    /* effeciencies are above the base running time on the OS,
		       the addition is due to multitasking needs,
		       real-time stuff */
		    bmipsGint = CPU.bmipsInt;
		    bmipsGFp = CPU.bmipsFp;
		    efficiencyG = 0.5; /* Operating system efficiency */
		    
		    Tmax =  85 degC; 
		    Tmin = -80 degC; // survival
		    
		    efficiencyD = 0.8; /* DSP specs */
		    bmipsDint = DSP.bmipsInt;
		    bmipsDfp = DSP.bmipsFp;
		  }
	      }
	      Component ConnectorBus{
		cost = 0.1;
		mass = 0.02*200 kg; /* about 4kg for 200 */
	      }
	      Component DataRecorder{
		cost = 0.8;
		mass = 2 kg;
		distrib power = distribution(max=15,avg=4,min=1);
	      }
	      Component Filters{ /* Analog filters & signal conditioners */
		cost = 0.1;
		mass = 1 kg;
	      }
	      Component SafetyCircuit{
		cost = 0.2;
		mass = (old rover.power/400 Watt) * 1 kg + 0.1 kg;
		volume = 10 cm * 10 cm * 2 cm;
	      }
	      Component Connectors
		{
		  Nconnectors = sum(key N of type sensor of Rover)*4;
		  mass = Nconnectors*(locomotion.L + locomotion.W)/(2 m)
		    * 50 gram; // ? is 50 gram/meter valid on average?
		  reliability = 0.99^Nconnectors;
		  /* this may be conservative as each "connector" may be
		     made of several strings - eg a ebox boundary */
		}
	      
	      mipsGfp = 50 * 10 ips;
	      mipsGint = 30 * 10 ips; /* 10 checks & adjustments per second */
	    } /* Computronics */
	  
	  Component Comm
	    {
	      
	      transparent mass = sum(key mass of Comm);
	      tiltMax = mission.slopemax + mission.tiltNormal 
		      + mission.latitude;
	      CommFreq = 8.9E9 Hz; 
	      Lambda = light.speed/CommFreq;
	      
	      CommGainSpace = -92.44 -20*log10(planet.earthdist/ (1 m)) 
		-20*log10(CommFreq/(1 Hz));
	      // The above is to be confirmed ?
	      CommGainTemp = 228.6 - 10*log10(mission.Tmax/(1 kelvin));
	      CommGainDataRate = -10*log10(downlink/(1 bps));
	      
	      Gain = sum(key "commgain*" of LunarRover) : 0.0;
	      
	      linkoverhead = 1.3; /* overhead on bandwidth for packetizing,
				     error correction, etc. */
	      transparent downlink = sum(key downlink of rover) * linkoverhead;
	      
	      CommGainEbNo = -2; // requirement for Nyquist QPSK, 2/3 turbo
	      CommGainMargin = -5; // Margin of 5 dB on top of all else
	      
	      Component Transmitter{} //  : AnnularLink {}   /* need to label EBox parts ? */
	      Component Receiver  : DSNDish {}
	      mipsGint = 100 * 10 ips;
	      /* 10 moniters per second to adjust bandwidth
		 and watch for blackouts */
	      mipsGfp = 100 * 10 ips;
	      
	      Gain >= 0;
	      
	      Component Messaging // is it also ? Packetizing & Error coding
		{
		  freq = 10 Hz;
		  mipsGint = freq * 100000 byte;
		}
	    } /* Comm(S. pole to Earth) */
	  
	  Component Sensors
	    {
	      type inherited sensor;
	      transparent mass = sum(key mass of sensors);

	      Component StarTracker /* Might need to make it redundant */
		{
		  type insulated;
		  resolution = 256 * 256 byte;
		  freq = 0.2 hz;
		  mass = 0.6 kg;
		  area = 5* 10 cm * 20 cm;
		  power = 3.25 watt;
		  cost = 0.75;
		  /* find three brightest in image, match in map */
		  mipsDint = resolution * 10 * freq; /* includes frame grab */
		  mipsGint = resolution * 2  * freq;
		  mipsGfp  = 100 kbyte * freq;
		}
	      Component IMU
		{
		  type inEBox;
		  readfreq = 400 Hz;
		  freq = 20 Hz;
		  mass = 0.7 kg;
		  volume = 10 cm * 10 cm * 10 cm;
		  power = 5 watt;
		  cost = 0.25;
		  mipsDfp = 8*8 * 10 byte * readfreq + 12*12 * 20 byte * freq;
		  /* read input at 400 Hz into a kalman filter,
		     use a larger filter at 20 Hz for output */
		}
	      
	      FovHneed = 360 degree; FovVneed = 269 degree; // Temporary ?
	      Component Imager  : CameraPanosphere{}
	      
	      Component SunTracker{ mass = 0.3 kg;}
	      Component TerrainSafeguarding
		{
		  type insulated;
		  
		  driveFOV = 120 degree; /* what do i need to drive */
		  SafetyFactor = 2.0;
		  
		  int Nfront = driveFOV/RadarMIR.fov*SafetyFactor;
		  /* *2 for double coverage */
		  Nside = (180 degree-driveFOV)/30 degree;
		  
		  Res = 2.0*min(wheel.d/2,wheel.w/2, mission.obstaclemax/2);
		  /* OK? takes care of high centering,
		     wheel trap in hole & climbable obstacles */
		  
		  /* assume single CPU mode, optimize if untenable */
		  transparent latency = sum(key latency of terrainsafeguarding) : 0;
		  speed = mission.speed;
		  /* reaction time is computing/comm time + stopping distance
		     from issue of command */
		  Lookahead1 = 0.9*(old latency * speed + 
				      speed*speed/planet.g/2.0/soil.mu_static);
		  /* Min lookahead, min distance for reaction, no reason
		     to look closer than that, Rho(reaction) - Al Kelly */
		  H = LookAhead1; /* sensor height hack, OK ? */
		  
		  /* angle wrt vertical */
		  LookAngle1 = atan((Lookahead1+SOffset)/H);
		  footprint  = tan(Lookangle1 + RadarMIR.FOV)*H - LookAhead1;
		  /* footprint is coverage of sensor on plat plane */
		  
		  /* lookahead so that each gridcell is seen at least four times,
		     four is a hack, is this correct ?*/
		  DeltaL = speed/RadarMIR.freq;
		  Lookahead2 = LookAhead1 + DeltaL + 
		    if(Res/DeltaL > 4.0, 0 m, DeltaL*(4 - Res/DeltaL));

		  SOffset = 0.0 m; 
		  /* where is sensor wrt front of machine, +ve => ahead */
		  LookAngle2 = atan((Lookahead2+SOffset)/H);
		  
		  
		  /* frontal area of sensor for MLI sizing */
		  area = (locomotion.L + Locomotion.W)*2*20 cm; /* major hack ? */
		  cost  = 2; /* includes R&D */
		  
		  /* Constraints */
		  footprint >= lookahead2 - lookahead1;
		  
		  Component RadarMIR
		    {
		      freq = 20 hz; /* Hz */
		      FOV = 30 degree; /* 30 deg */
		      /* no of gridcells covered per sensor read,
		         compute area of "sensor triangle" upto lookahead2
			 and then divide by gridcell size */
		      cellscovered = lookahead1 * 
			((int) ((lookahead2 - lookahead1)/Res + 1.0))*
			tan(FOV)/2/Res;
		      /* for single coverage on sides */
		      N = Nfront + Nside*3; /* *3 to cover sides & back */
		      /* evidence grid maps computing = 10 instructions * no of
			 grid cells covered by a sensor * no of sensors *
			 sense freq - Hans Moravec */
		      mipsGint = 10 byte * SafetyFactor * freq * cellscovered;
		      latency  = old mipsGint/(board.cpu.bmipsint*boards.OSgOverhead) * 1 second
			+ 1.0/freq; /* one cycles to process & transmit to CPU mem*/
		      mass  = 0.4 kg * N;
		      power = 0.8 watt * N;
		      sensorangle = (LookAngle1+LookAngle2)/2.0;
		    }
		  Component Planning
		    {
		      MapRes = Res;
		      freq = speed/Res; /* is there something better? */
		      Narcs = 9;
		      mipsGfp = Narcs * (lookahead2-lookahead1)*1.5/MapRes
			* 50 byte * freq;
		      /* correct? Al Kelly */
		      
		      /* latency is processing time + delay due to freq diff
			 between this module & Radar,an extra 1/radarMIR.freq is
			 implicit due to possible bad synchronization */
		      latency = mipsGfp/board.cpu.bmipsfp/boards.OSgOverhead
			* 1 second
			+ if( (freq == radarMIR.freq),
			     0, (freq/radarMIR.freq)/radarMIR.freq);
		    }
		} /* TerrainPerception */
	      Component internalSensors
		{
		  type sensor, inEBox;
		  N = rover.Nactuators + 2; /* .5 for electronics box */
		  mass = 0.25 kg*N ;
		  power = 0.25 watt*N;
		  cost = 0.1 megadollar + 0.05 megadollar*N;
		}
	    } /* sensors */
	  
	  Component Software
	    {
	      Component SensorRead
		{
		  mipsGint = rover.Nsensors* 50 byte * 100 Hz; /* is this valid ? */
		}
	      Component telemetry  : MinimalTelemetry{}
	    }
	  
	  Component SciencePackage : SciencePackageNone{}
	  
	} /* Rover */
      
      /* Rename Lander into Mars98Lunar and push to options
	 Replace by Component Lander : Mars98Lander;
	 Component PhobosLander // Hack, tries to fit in two rovers
	 {
	 mass = 6100 kg; // including ramps
	 NRovers = 2;
         rover.mass*NRovers <= 500 kg; // ramps included in lander
	 H = 1.0m; // Verify ? this and CG
	 CGz = H/2.0; 
	 RoverCGL = Rover.L/2.0 - Rover.CGLoff;
	 RoverCGL <= 1.0 m; // OK?
	 LanderRoverMass = (old Rover.mass)*NRovers + mass;
	 LanderRoverCGz = ((RoverCGL+H)*(old Rover.mass)*Nrovers +
	    mass*CGz)/LanderRoverMass; // Rovers vertical on Phobos
	 cost = 40 megadollar;
	 }
     */
	 
      Component Lander
	{
	  transparent mass = sum( key mass of lander ) : 200 kg;
	  H = 0.4 m; /* Highly optimized lander, 0.0 if lander rover */
	  CGz = H/2.0;
	  LanderRoverMass = (old Rover.mass) + mass;
	  LanderRoverCGz = ((Rover.CGz+H)*(old Rover.mass) +
	    mass*CGz)/LanderRoverMass;

	  Component Rocket : RocketAnalysis {}
	  cost = 15;
	  Component Fuel
	    {
	      fuel_spent = old fuel_spent + fuel_rate( old total_thrust ) * delta_time :
	      0 kg, clock coriolis;
	      
	      fuel_rate( thrust ) = thrust / isp / 9.81 [m/s2] : 0 [kg/s], clock coriolis;
	      
	      Isp = 220 sec;
	      dVmidcourse = 50 [m/s]; /* m/s, for a 1.5 deg max error, confirm */
	      Mmidcourse = old spacecraft.mass - old dMmidcourse : 1500 kg;
	      dMmidcourse = rocket.DeltaM(Mmidcourse, Isp, dVmidcourse) : 35 kg;
	      
	      // dVacs = ; /* m/s, attitude control during SRM burn */
	      Macs = Mmidcourse - old dMacs;
	      dMacs = Mmidcourse/1456*25 : 25 kg; /* kg */ /* LOKI -- needs ival */
	      
	      dVterminal = 200 [m/s]; /* m/s */
	      Mterminal = Macs - old SRM.mass - old dMterminal;
	      dMterminal = rocket.DeltaM(Mterminal, Isp, dVterminal) : 50 kg;
	      
	      Mhover = Mterminal + dMterminal;
	      hangtime = 90 s; /* seconds of course */
	      /* ratio of fuel mass/s to craft mass */
	      dMhover = rocket.DeltaMhover(Mhover,Isp,hangtime,planet.g) : 15 kg;
	      /* kg, for hovering above surface */
	      
	      /* dMmargin = 0.0*mass : 25 kg;  kg, fuel margin */
	      
	      mass = sum (key "dm*" of fuel) - fuel_spent : 130 kg, 
	      clock coriolis, main_clock;
	    }
	  Component Propulsion
	    {
	      Component DE 
		{
		  mass = 60 kg; /* Includes ACS & tubing as of now */
		} /* Descent Engine */
	      Component ACS {} /* Attitude Control System */
	      Component BellyCamera 
		{
		  mass = 2 kg;
		  cost = 0.25;
		}
	      Component Radar /* landing radar */
		{
		  mass = 3 kg;
		  cost = 0.5;
		}
	      Component Software {}
	      Component Disconnect {}
	      Component Structure {}
	      Component Electronics {}
	    } /* propulsion */
	  
	  shape body_shape = shape "booster_shape.txt" * ( 4.8, 1, 4.8 );
	  triple position = triple( 0, 12.5, 0 ) kilometer;
	  triple velocity = -triple( 0, 50, 0 ) [m/s];
	  altitude = dot( y_axis, position ) - 1 m, clock coriolis;
	  
	  weight = ( mass - srm.mass ) * planet.g, clock coriolis;
	  
	  /* gravity */
	  link lander : force weight * -y_axis, triple( 0, 0, 0 );
	  
	  Component Thrust_Control/*, clock coriolis*/
	    {
	      /* target velocity for spacecraft */
	      real vtarget = if( altitude > 2 kilometer, -150 [m/s],
				altitude > 1 kilometer,
				-( ( altitude - 1 kilometer ) * 0.1 / 1 s
				  + 50 [m/s] ), altitude > 5 m,
				-( altitude - 5 m ) / 20 s - .25 [m/s],
				-0.25 [m/s] ), clock coriolis;
	      
	      /* if dropping too quickly then delta_v negative */
	      real delta_v = dot( velocity, y_axis ) - vtarget,
	      clock coriolis;
	      
	      /* dv1 >= 0 */
	      dv1 = if( delta_v >= 0 [m/s], 0 [m/s], -delta_v * 1.2 ), clock coriolis;
	      /* error < 1 */
	      error = if (delta_v >= 0 [m/s], 0 [m/s], old error * .7 + dv1 * .1) : 0 [m/s],
	      clock coriolis;
	      
	      /* desired_thrust >= 0? */
	      thrust = ( dv1 * 1.2 [s/m] + error * 1.4 [s/m] ) * weight,
	      clock coriolis;
	      
	      transparent maxThrust = 3300 newton, clock coriolis;
	      
	      /* real thrust values in newtons */
	      w_thrust = max( min( maxThrust, thrust / 6 ), 0 newton ), clock coriolis;
	      e_thrust = max( min( maxThrust, thrust / 6 ), 0 newton ), clock coriolis;
	      nw_thrust = max( min( maxThrust, thrust / 6 ), 0 newton ), clock coriolis;
	      ne_thrust = max( min( maxThrust, thrust / 6 ), 0 newton ), clock coriolis;
	      sw_thrust = max( min( maxThrust, thrust / 6 ), 0 newton ), clock coriolis;
	      se_thrust = max( min( maxThrust, thrust / 6 ), 0 newton ), clock coriolis;
	      total_thrust = sum( key "*_thrust" of Thrust_Control ) : 0 newton,
	      clock coriolis;
	      
	      /* thrusters -- two sets of three.  center thruster lies in the xy plane,
		 at a 10 degree angle from the negative y axis.the side two are 10
		 degrees from the negative y axis in the xy plane and 15
		 degrees from the negative y axis in the yz plane */
	      x_thruster_angle = 10 degree;
	      z_thruster_angle = 15 degree;
	      
	      x_tan = tan( x_thruster_angle );
	      z_tan = tan( z_thruster_angle );
	      
	      y_side = 1 / sqrt( x_tan * x_tan + 1 + z_tan * z_tan );
	      x_side = x_tan * y_side;
	      z_side = z_tan * y_side;
	      
	      y_main = 1 / sqrt( x_tan * x_tan + 1 );
	      x_main = x_tan * y_main;
	      
	      /* thrusters are labeled east_wast if in the pos/neg * direction, or
		 north_south for the pos/neg z direction */
	      triple thrust_west = gvec( x_main, y_main, 0 ), clock coriolis;
	      triple thrust_east = gvec(-x_main, y_main, 0 ), clock coriolis;
	      
	      triple thrust_northwest = gvec( x_side, y_side, -z_side ), clock coriolis;
	      triple thrust_northeast = gvec(-x_side, y_side, -z_side ), clock coriolis;
	      triple thrust_southwest = gvec( x_side, y_side, z_side ), clock coriolis;
	      triple thrust_southeast = gvec(-x_side, y_side, z_side ), clock coriolis;
	      
	      triple wthrust_bp = triple( -2.4, -.30, 0 ), clock coriolis;
	      /* -.50, -.30 scaled */
	      triple ethrust_bp = triple( 2.4, -.30, 0 ), clock coriolis;
	      /* .50, -.30 scaled */
	      
	      /* forces due to thrusters */
	      link spacecraft : force w_thrust * thrust_west, wthrust_bp;
	      link spacecraft : force e_thrust * thrust_east, ethrust_bp;
	      link spacecraft : force nw_thrust * thrust_northwest, wthrust_bp;
	      link spacecraft : force ne_thrust * thrust_northeast, ethrust_bp;
	      link spacecraft : force sw_thrust * thrust_southwest, wthrust_bp;
	      link spacecraft : force se_thrust * thrust_southeast, ethrust_bp;
	    } /* thrust control */

	  Component LanderRoverinterface {mass = 0.005*(old Rover.mass);
					  cost = 0.2;}
	} /* Lander */


      /* 
	 Rename the SRM as ThiokolStar & push into options
	 Rplace by Component SRM : ThiokolStar;
	 
	 Add option
	 Component SRMnone { mass = 0.0kg; CGz = 0.0m; H= 0.0m;};
      */
            Component SRM
	{ /* Solid Rocket Motor */
	  transparent mass = sum(key mass of SRM) : 974 kg;
	  CGz = H/2.0;
	  H = 1.1 m;
	  Cost = 250 [dollar/kg] *old mass*2.2;
	  Component Fuel
	    {
	      dV = 2560 [m/s]; /* m/s */
	      /* mass = rocket.deltaM(old spacecraft.Mass,
		 290 s, dV); */
	      mass = 861*old spacecraft.mass/ 1500;
	    }
	  Component drySRM /* Thiokol Star 37FXP */
	    {
	      mass = 63.6 kg;
	    }
	  Component LanderInterface
	    {
	      // Is the magic number 0.002 or 0.005
	      Component Clamp {mass = 0.002*(old Rover.mass + lander.mass);}
	      Component Harness {mass = 2 kg;} 
	      cost = 0.2;
	    }
	} /* SRM */

    } /* Spacecraft */
  
  Component DSN {cost = 1500 [dollar/hour] * duration;}
  
  Component SystemEngineeringDDT
    {
      overhead = 2.0;
      Component Coreteam
	{
	  personCost = 80000 dollar *Overhead;
	  cost = personCost*work; /* 80k / engineer */
	  work = 20*3;
	}
      Component facility {cost = 6; }
      Component IronBird {cost = 2; }
      Component prototype {cost = 2; }
      Component mobilityTestbed {cost = 0.6; }
      Component testing
	{
	  Component equipment {cost = 1; }
	  Component facility {cost = 0.5; }
	  Component radiationFacility {cost = 0.2; }
	  Component transportation {cost = 0.6; }
	  Component burnInrisk {cost = 3; }
	} /* testing */
      
      Component Management 
	{
	  cost = 0.15*sum(key cost of LunarRover 
			  except set( management.cost, contingency.cost ) );
	}
      Component Contingency {cost = 0.2*management.cost/0.15*1.15;}
      
      Component Software 
	{
	  unitcost = 65000 dollar * overhead;
	  cost = unitcost*sum(key work of LunarRover); /* 60-70 + overhead */
	  Component Simulation {work = 5*2;}
	  Component CAD {work = 15*2;}
	  Component Ground {work = 2*2;}
	  Component Science {work = 2*2;}
	  Component Test {work = 5*2;}
	  Component Injection {work = 5*2;}
	  Component Satellite {work = 2*2;}
	  Component Lander {work = 5*2;}
	  Component HousekeepingSafeguarding {work = 5*2;}
	  Component Contingency {work = 5*2;}
	} /* software */
      
    } /* systemEngineeringDDT */
  
} /* LunarRover */

/********************** Optional (sub) systems & components ****************/


Component Options
{
  dormant = true;
  Component Parts
    {
      Component Comm
	{
	  Component PhasedArray
	    {
	      HPBW = [cast degree] 58.5*comm.lambda/D;
	      factor = 1./(1. + sin(comm.tiltMax + 1.5*HPBW)) : 1;
	      dx = old factor*comm.lambda;
	      dy = old factor*comm.lambda/sqrt(3.); // why is dy different than dx
	      Nelements = min((old Nelements)*10^((-(old comm.gain)/4.0 + 0.1)/10), 750):683;
	      D = 2.*sqrt(Nelements*2.*dx*dy/Pi);
	      
	      PowerElement = min(old Powerelement * 
				 (Nelements/683.0), 25 milliwatt)
	      : 18.8 milliwatt;

	      PowerTransmit = Nelements*PowerElement;
	      Pefficiency = 0.126;
	      Power = PowerTransmit/Pefficiency + 6.0 watt;
	      /* ~100W if 638 elements at 18.8 mW, 0.126 efficiency */
	      Mass = 12.8 kg/683*Nelements + 2.0  kg; /* Check ? */
	      efficiency = 0.7;
	      CommGainArray = 10*log10(Nelements);
	      CommGainLoss = -10*log10(efficiency) + 3.0;
	      CommGainPow = 10*log10(PowerTransmit/1 watt);  // ?
	      CommGainLineloss = -3.0; // Line Loss, Gl
	      CommGainPointing = -10*log10(cos(tiltmax) );

	      PowerElement <= 25 milliwatt;
	      Nelements <= 683;
	      
	      Component mast  : EmptyMast{}
	    }
	  
	  Component OmniLink
	    {
	      mass = 2.1 kg;
	      distrib power = distribution(avg = 0 , min = 0, max = 40) watt;
	      Component mast  : EmptyMast{}
	    }
	  
	  Component ParabolicAntenna
	    {
	      CommGainDish = 20*log10(Pi*D/Comm.Lambda) : 28; // 0 dB => omni
	      CommGainDishLoss = -10*log10(efficiency) + 3.0;
	      Component mast  : EmptyMast{}
	      
	      D = 10^((-Comm.gain + old CommGainDish + 0.1)/20)*Comm.Lambda/Pi;
	      HPBW = [cast degree] 58.5*comm.lambda/D;
	      
	      PowerTransmit = 15 Watt;
	      Pefficiency = 15 Watt / 40 Watt;
	      Power = PowerTransmit/Pefficiency + 5.0;
	      /* ~100W if 638 elements at 18.8 mW, 0.126 efficiency */
	      MassDish = 1 kg * (D^2/(0.5 m * 0.5 m));
	      Mass = 4 kg + MassDish;
	      efficiency = 0.55;
	      CommGainPow = 10*log10(PowerTransmit/1 watt);  // ?
	      CommGainLineloss = -3.0; // Line Loss, Gl
	      
	      PowerTransmit <= 15 Watt;
	      
	      Component Pointing
		{
		  mass = if(HPBW*1.5 < comm.tiltmax, 
			    3 kg * (MassDish/1 kg),0.0 kg);
		  power = if(HPBW*1.5 < comm.tiltmax, 
			     20 Watt * (MassDish/1 kg), 0.0 kg); 
		  // Deepak?
		}
	    } /* ParabolicAntenna */
	  
	  Component AnnularLink
	    {
	      /* For omni comm from South Pole */
	      /* If there is no satellite */
	      Component transciever
		{
		  HPBW = 360 degree;
		  mass = 5 kg;
		  distrib power = distribution(max=40,avg=40,min=0) watt;
		  efficiency = 0.3;
		  TransmitPower = power * efficiency;
		  CommGainPower = 10*log10(TransmitPower/(1 Watt)); // ?
		  CommGainAnt = -2; // For about 50 deg annular
		}
	      Component Antenna
		{
		  type mastpayload;
		  mass = 1 kg;
		  CommGain = 8; // ?
		}
	      Component Cable
		{
		  type mastpayload;
		  H = mast.h + 1.0 m; /* meters */
		  mass = H*0.5 [kg/m];
		  CommGain = -0.2*H/(1.0 m); // ?
		}
	      Component backup
		{
		  mass = 3 kg;
		  distrib power = distribution(max=30,avg=0,min=0);
		}
	      /* This Comm mast, mass is a function of height and mast payload */
	      Component Mast{
		H = 6 m;
		Component payload{transparent mass = sum(key mass of type
							 mastpayload of rover);}
		mass = payload.mass/5.0*H/(3.0 m)*(Earth.g/planet.g);
		/* add 0.2 mass every 3 m of "earth height" ? */
	      }

	      Component Pointing
		{type mastpayload;

		  mass = if(HPBW*1.5 < comm.tiltmax, 3 kg * D/0.5 m,0.0 kg);
		  power = if(HPBW*1.5 < comm.tiltmax, 20 Watt * D/0.5 m, 0.0 watt); // Deepak?
		}
	    } /* AnnularLink */
	  
	  Component EmptyMast{H = 0 m; mass = 0 kg;}
	  
	  Component DSNDish
	    {
	      D = 34 m; // 22 m and 34 m are available

	      CommGainWaveDistort = -0.2;
	      CommGainModulate = -0.2;
	      CommGainPolarise = -0.2;
	      CommGainPoint = -0.2; // Hot Damn
	      CommGainCircuit = -0.5;
	      CommGainPhase = -0.2;
	      CommGainFeed = -0.1;
	      efficiency = 0.55; /* Should DSN not be better 0.6-0.7? */
	      CommGainDish = 20*log10(Pi*D/Comm.Lambda) + 10*log10(efficiency) + 3.0; // 3.0 OK?
	      
	    }
	  
	}
      
      Component Power
	{
	  Component SiliconCell
	    {
	      unitcost = 800 [dollar/watt];
	      efficiency = 0.14 - 0.02; /* 0.02 launch and attrition */
	    }
	  Component GaAsCell
	    {
	      unitcost   = 1200 [dollar/watt];
	      efficiency = 0.18 - 0.02;
	    }
	  
	  Component RTG // Note, there is an additional 20% structural margin
	    {
	      powerOut = 1.1* powerNeed;
	      nightpowerOut = powerOut;
	      Component AMTEC
		{
		  mass = 1.0 [kg/watt]/40 * powerOut;
		  efficiency = 0.15;
		  Tcold = 550 Kelvin; /* max cold temp AMTEC can reach */
		  Thot = 750 Kelvin;
		}
	      Component GPHS
		{
		  output = powerOut/AMTEC.efficiency;
		  heatflow = output*(1.0 - efficiency);
		  mass = output/(250 watt) * 1.5 kg; /* use ceil function */
		}
	      Component Thermal /* include thermal ?, MLI, heatpipe, radiator */
		{
		  mass = 19.0/33.0 * (AMTEC.mass + GPHS.mass);
		  // To give about 60 kg for a 450 W system
		}
	    }
	  
	  Component SolarArray
	    { /* 4 mil cell, 3 mil coating, single sided */
	      Component Cell  : SiliconCell{}
	      W = 2.5 m;
	      H = 2.0 * powerNeed/ insolation/ efficiency/W;
	      /* is 2.0 OK? */
	      area = H * W; /* w x h */
	      mass = 0.75 [kg/m2] * area * 2; /* 2 - structural overhead */
	      powerOut(itheta = 0) = if( (cos(itheta) < cos(45 deg)), 0.0 watt,
					area* insolation * efficiency * cos(itheta) ); /* units LOKI */
	      nightpowerout = 0.0 Watt;
	      cost = powerOut*unitcost + 0.5 megadollar; /* .5 for DD&T */
	      powerOut >= old Rover.power;
	      
	      /* Bug ? default type ips is not attached if the ips is omitted */
	      mipsGfp = 100 * 10 ips; /* need to do constant moniter and regulate */
	      
	      Component Actuator
		{ /* Solar array 1 axis pointing */
		  mass = 4 kg;
		  distrib power = distribution(max=15,avg=5,min=0);
		  cost = 0.5 megadollar;
		}
	    }
	  
	  Component Battery /* check data against real numbers ? */
	    /* incorporate life (total usage time), power density, discharge/charge rates
	       Since packaging matters energy desnity is also function of total
	       charge, where as size incr, density increases */
	    {
	      Component AgZnHR190
		{ /* HR190 - Yardney */
		  unitMass = ( 1 [kg/whr] ) / 190; /* kg/Whr */
		  unitVolume = 0.256834 [cc/Whr]/unitmass; /* cc/Whr ?*/
		  Tmax = 40 degC;
		  Tmin = -20 degC;
		  
		  d1 = 0.2; n1 = 15*n2; /* refer to Cataldo's e-mail in LRI/discovery 200, 242, ... */
		  d2 = 0.9; n2 = 30; /* these numbers are so subjective 25-50 ? */
		  Life(e) = 10.0^(log10(n2/n1)*log10(e/energyNeed/d1)/ log10(d2/d1)
				  + log10(n1));
		  
		  mass = energyNeed*unitMass;
		  volume = mass*unitVolume; /* cc */ /* CHANGED energy to energyNeed */
		  cost = 0.5 megadollar + 1.5 megadollar*mass/40 kg; /* includes engineering, DD&T, proto */
		  efficiency = 0.8; /* charging effectiveness */
		}
	      
	      Component SuperNiCd
		{
		  unitMass = 1 [kg/Whr]/60.0;
		  unitVolume = 1.2 liter/kg;
		  Tmax =  50 degC;
		  Tmin = -40 degC;
		  
		  mass = energyNeed*unitMass;
		  volume = mass*unitVolume; /* cc */ /* CHANED energy to energyNeed */
		  d1 = 0.2; n1 = 40000;
		  d2 = 0.9; n2 = 3000;
		  Life(e = 0 joule) = 10.0^(log10(n2/n1)*log10(e/energyNeed/d1)/ log10(d2/d1)
				  + log10(n1));
		  cost = 0.5 megadollar + 1.5 megadollar*mass/(40 kg); 
		  /* ? includes engineering, DD&T, proto */
		  efficiency = 0.8; /* charging effectiveness */
		}
	      
	      Component NiH2 // Nickel Hydrides
		{
		  unitMass = 1 [kg/whr]/70.0;
		  unitVolume = 1.2 [liter/kg];
		  Tmax =  10 degC;
		  Tmin =   0 degC;
		  
		  mass = energyNeed*unitMass;
		  volume = mass*unitVolume; /* cc */ /* CHANED energy to energyNeed */
		  d1 = 0.2; n1 = 70000;
		  d2 = 0.9; n2 = 8000;
		  Life(e) = 10.0^(log10(n2/n1)*log10(e/energyNeed/d1)/ log10(d2/d1)
				  + log10(n1));
		  cost = 0.5 megadollar + 1.5 megadollar*mass/(40 kg); /* ? includes engineering, DD&T, proto */
		  efficiency = 0.8; /* charging effectiveness */
		}
	      
      	      Component LiIon
		{
		  unitMass = 1 [kg/whr]/150.0;
		  unitVolume = 1.2 [liter/kg];
		  Tmax =  50 degC;
		  Tmin = -10 degC;
		  
		  mass = energyNeed*unitMass;
		  volume = mass*unitVolume; /* cc */ /* CHANED energy to energyNeed */
		  d1 = 0.2; n1 = 40000;
		  d2 = 0.9; n2 = 3000;
		  Life(e) = 10.0^(log10(n2/n1)*log10(e/energyNeed/d1)/ log10(d2/d1)
				  + log10(n1));
		  cost = 0.5 megadollar + 1.5 megadollar*mass/(40 kg); 
		  efficiency = 0.8; /* charging effectiveness */
		}
	      
	    } /* Battery */
	}
      
      Component Thermal
	{
	  Component MLI     /* N & A need to be defined */
	    {
	      sigma = 5.67E-8 [kg/s3/kelvin4]; 
	      MLImass(area,layers) = area*layers*(0.04 [kg/m2] + 0.01 [kg/m2]);
	      /* 1 mil (25 um - 1600 kg/m^3) kapton, 0.01 extra for
		 stiching n such */
	      mass = MLImass(area,nlayers);
	      E = 0.05/Nlayers*2; /* 2 is factor of safety */
	      area = 5.0 m2; nlayers = 10;
	      As = Area/3.0; /* assuming cubic shape, 1/4 if spherical */
	      heatflow(Ti,To, sarea=area) = sigma*SArea*E*(Ti^4-To^4); /* outflow */
	    }
	  
	  Component Radiator  /* Specify Area */
	    {
	      Component surface  : SilverTeflon{}
	      
	      sigma = 5.67E-8 [kg/s3/kelvin4];
	      /* L = (Tr > 310 kelvin)?0.9:0.1 louvre position 0.9 to 0.1 */
	      heatflow(Ti,To,itheta=0) =
		sigma*Area*(if( (Ti > 310 kelvin),0.9,0.1 ))*surface.Eir*(Ti^4-To^4) -
		  insolation*Area*surface.Es*cos(itheta); /* outflow */
	      SArea(Ti,To, heat, maxi) = heat/(sigma*(if((Ti > 310 kelvin),0.9,0.1))*
					 surface.Eir*(Ti^4-To^4) -
					 insolation*cos(maxi)*surface.Es); /* Works for worst case sun angle */
	      HeatCold(Ti,To, sarea=area) =
		sigma*sArea*(if((Ti > 310 kelvin),0.9,0.1))*surface.Eir*(Ti^4-To^4);
	      mass = Area*2.2 [kg/m2]; /* real values?, fn of heat flow? */
	      cost = surface.costDDT + surface.unitcost*Area;
	    }
	  
	  Component SilverTeflon /* ? */
	    {
	      Eir = 0.85;
	      Es  = 0.25;
	      costDDT  = 0.1 megadollar;
	      unitcost = 0.4 [megadollar/m2];
	    }
	  
	  Component Black /* ? */
	    {
	      Eir = 0.85;
	      Es  = 0.85;
	      costDDT  = 0.1 megadollar;
	      unitcost = 0.15 [megadollar/m2];
	    }
	  
	  Component ConvectiveHeatPipe
	    {
	      dhdt = 100 [watt/kelvin];
	      dT(heatflow = 0) = heatflow/dhdt;
	      D = 0.1 m; /* dia is probably less, give it room */
	      mass = H/(2 [m/kg]); /* 1 kg/ 2 m ?*/
	      volume = Pi/4.0*H*D*D;
	      area = H*D*Pi;
	    }
	  
	  Component Heaters{
	    type inEBox;
	    mass = 0.5 kg; /* make this a fn of heat output */
	    volume = 5 cm * 5 cm * 5 cm;
	  }
	  
	  Component HeatSourceRHU 
	    {
	      unitmassRHU = 60 gram/ 2 Watt;
	      unitmassGPHS = 1.5 kg/ 250 Watt;
	      
	  
	      /*? Lalit's Hack to avoid dependency problem, Original Code: */
	      // nightPowerOut = if(mission.night > 0.0 sec,
			//	 nightPowerNeed - sum(key nightPowerOut of rover), 0.0 watt);
	      /* Now the Hack which should be replaced by above code once
		 the bug is fixed, damn bugs ? */

	      Hack = nightPowerNeed - PowerSource.nightpowerOut;
	      nightPowerOut = if(mission.night > 0.0 sec,
				 Hack, 0.0 watt);
	      /* END hack */

	      mass = min((real) ((int) (nightPowerOut * unitmassGPHS + 0.5 kg)),
			 nightPowerOut * unitmassRHU);
	    }
	  
	  Component HeatSourceNone { nightPowerOut = 0.0 Watt; }
	  
	}
      
      Component Propulsion
	{
	  Component Launcher
	    {
	      Component DeltaII7925
		{
		  timelead = 3 year;
		  // Actual formula, BUG confirm that factor and factor2 are the same ?
		  massTLI = 1315 kg * (1.0 - 1.4/2.3 * if(spacecraft.CGz > 0.5 m,
							  (spacecraft.CGz - 50 cm)/1.143 m, 0.0)); 
		  factor = if(spacecraft.CGz > 0.5 m,
			      (spacecraft.CGz - 50 cm)/1.143 m, 0.0);
		  factor2 = if(spacecraft.CGz > 50 cm,
			       (spacecraft.CGz - 50 cm)/1.143 m, 0.0);

		  Dia = 2.3 m; /* 2.3 meter payload fairing */
		  H = 1.8 m;
		  dia2 = 0.0 m;
		  H2 = 3.0 m;
		  InjectionError = 1.5 degree; /* deg */
		}
	      Component Ariane5miniSylda
		{
		  timelead = 4 year;
		  distribution massTLI = distrib(min=850,max=1200,avg=1200) kg;
		  Dia = 3.6 m; /* 2.3 meter payload fairing */
		  H = 1.8 m;
		  dia2 = 0.0 m;
		  H2 = 0.0 m;
		  InjectionError = 2.0 degree; /* confirm? */
		}
	    } /* launcher */
	  
	  Component RocketAnalysis
	    {
	      G = 9.81 [m/s2]; /* m/s^2 */
	      DeltaM(Mf, Isp, dV) = Mf*(exp(dV/(Isp*earth.G))-1)/exp(dV/(Isp*earth.G));/* kg */
	      /* confirm? */
	      /* Returns fuel needed given final Mass after burn, Isp & deltaV
		 Do not forget to include rocket dry mass in your Mf, may need
		 iterations, so what does not?
		 
		 dV: Delta V
		 Vj: Jet velocity, Vj = Isp*G;
		 Mf: Spacecraft mass
		 mu = exp(dV/Vj)
		 Propellant required (approx) = dV/Vj*mo
		 Propellent required (actual) = (mu - 1)/mu * mf
		 Isp: lbf/(lbmass/sec) [in Earth pounds force]
		 */
	      
	      DeltaMhover(Mf, Isp, hangtime, grav) = Mf/(1/(grav*hangtime/Isp/Earth.g) - 1.0);
	      
	    } /* Rocket */
	  Component ThiokolStar37XFP
	    {
	      massdry = 73 kg;
	      massfuel = 884 kg;
	      mass = massfuel + massdry;
	      life = 66.5 s;
	      thrust = 2097 kiloNewton/life;
	      dia = 1.51 m;
	      H = .93 m; /* confirm? */
	    }
	} /* Propulsion */
      
      Component Mechanism
	{
	  Component RigidWheel  /* need to specify n no of wheels,d dia,w width*/
	    { /* rigid wheel model */
	      
	      n = no; /* n wheels */
	      
	      d = old dia; /* m */
	      w = old wid; /* m */
	      
	      Contact(load) = sqrt(d^2/4 - (d/2 - sinkage(wgt=load) )^2 );
	      contactArea(load) = Contact(load)*wheel.w;
	      pressureMax = (old Rover.mass)*planet.G/
		(2*contactArea(WheelLoadMax));
	      /* load is wheel load for given robot weight and slope */
	      load(theta = 0,wgt= Rover.weight) = wgt/n*cos(theta)*1.1;
	      /* 1.1 accounts for imbalance, cross slopes */
	      
	      SinkMax = sinkage(0, load(0 radian, locomotion.weight));
	      /* for rigid wheel, sinkage equation */
	      /* Bad hacking to get around bugs in unitcasting & unit checks
		 Fix at some point ? */
	      sinkage(theta,wgt=locomotion.weight) = 
		((3*load(theta,wgt)/(1 Newton)/((3-soil.n)* 
				 (soil.kc+ W * soil.kphi)/ (1 [Newton/m2])
				      * sqrt(D/(1 m))))^(2/(2*soil.n+1))) *
					(1 meter);
	      
	      /* Soil thrust */
	      H(theta,wgt=locomotion.weight) = 
		(soil.c*contactArea(load(theta,wgt)) + 
		 load(theta,wgt)*tan(soil.phi));
	      
	      
	      /* Compaction resistance */ /* LOKI fractional power of D */
	      /* 	  Rc(theta,wgt) = ((3*load(theta,wgt))^((2*soil.n+2)/(2*soil.n+1)) */
	      /* 		       /sqrt(D^((2*soil.n+2)/(2*soil.n+1)))) */
	      /* 	    / ((3-soil.n)^((2*soil.n+2)/(2*soil.n+1))*(soil.n+1) */
	      /* 	       *(soil.kc+W*soil.kphi)^(1/(2*soil.n+1))); */
	      Rc(theta,wgt) = (1 Newton) * 
		(3*load(theta,wgt)/(1 Newton)/sqrt(D/(1 m)))^((2*soil.n+2)/(2*soil.n+1))
		/ ((3-soil.n)^((2*soil.n+2)/(2*soil.n+1))*(soil.n+1)
		   *((soil.kc+W*soil.kphi)/(1 [Newton/m2]))^(1/(2*soil.n+1)));
	      
	      /* Buldozing resistance */
	      alpha = 1.0 [Newton/m3];

	      /* Unitcasting hacked here too ? */
	      Rb(theta,wgt) = (0.5*alpha*W*sinkage(0,load(theta,wgt))^2*
			   tan(pi/4 + soil.phi/2)^2)
		+ 2*soil.c*W*sinkage(0,load(theta,wgt))*
		  tan(pi/4+soil.phi/2);
	      
	      /* Rolling resistance */
	      Rr(theta,wgt) = soil.fr*load(theta,wgt);
	      
	      /* Gravitational resistance */
	      Rg(theta,wgt=locomotion.weight) = wgt/n * sin(theta);
	      
	      /* Acceleration resistance */
	      Ra(a=0 [m/s2]) = (old Rover.mass)/n * a;
	      
	      /* Obstacle climbing resistance */
	      Ro(size=0) = 0;
	      
	      /* Total resistance excluding obstacle, steer, acceleration */
	      Rall(theta=0 radian,wgt=locomotion.weight) = Rc(theta,wgt)+
		Rb(theta,wgt) + Rr(theta,wgt) + Rg(theta,wgt);
	      
	      DP(theta=0 radian,wgt=locomotion.weight) = 
		H(theta,wgt) - Rall(theta,wgt);

	      /* Driving Torque */
	      
	      defl = 0.0 m; /* deflection of wheel */
	      Tw(slope=0 radian,wgt=rover.weight) = Rall(slope,wgt)*(D/2-defl);
	      

	      Torquemax = Tw(mission.slopemax, wgt = rover.weight);
// replaceafter deBUG? wgt=rover.weight); /* LOKI - NAN */ - seems fixed whee...
	      
	      NCycles = mission.distance/(Pi*D);
	      
	      /* Wheel Design Constraints */
	      sinkage(theta=0 radian,wgt=WheelLoadMax) <= 0.1*wheel.D;
	      constraint sinking = pressureMax <= soil.pressureMax;
	    } /* Wheel */
	} /* mechanism */
      
      
      /* Need to add rad tolerance knowledge in */
      Component Computers /* perhaps Spec ratings will also be useful ? */
	{
	  component mipsR10000
	    {
	      type computer;
	      bmipsInt = 428 mips;
	      bmipsFp = 760 mips;
	      bpower = 6 watt; /* ? */
	      dcost  = 0.75 megadollar; // development cost
	      Tmax =  70 degC;
	      Tmin = -80 degC;
	    }
	  component mipsR3000
	    {
	      type computer;
	      bmipsInt = 28 mips;
	      bmipsFp = 36 mips;
	      bpower = 6 watt;
	      dcost  = 0.1;
	      Tmax =  85 degC;
	      Tmin = -80 degC;
	    }
	  component PowerPC750
	    {
	      type computer;
	      bmipsInt = 496 mips;
	      bmipsFp  = 336 mips;
	      bpower   = 5.7 watt;
	      dcost  = 0.75;
	      Tmax =  70 degC;
	      Tmin = -80 degC;
	    }
	  component PentiumPro
	    {
	      type computer;
	      bmipsInt = 292 watt;
	      bmipsFp = 248 watt;
	      bpower = 23 watt;
	      dcost  = 1.25;
	      Tmax =  60 degC;
	      Tmin = -80 degC;
	    }
	  component TIC60
	    {
	      type computer;
	      bmipsInt = 1600 mips;
	      bmipsFp  = 1600 mips;
	      bpower   = 11 watt; /* ? */
	      dcost    = 0.75 megadollar+ 0.5 megadollar; /* software lib dev costs */
	      Tmax =  70 degC;
	      Tmin = -80 degC;
	    }
	  component TIC80
	    {
	      type computer;
	      bmipsInt = 2200 mips;
	      bmipsFp  = 2200 mips;
	      bpower   = 11 watt;
	      dcost  = 0.75 megadollar + 0.5 megadollar; /* 0.5 - software lib dev costs */
	      Tmax =  70 degC;
	      Tmin = -80 degC;
	    }
	  component Empty
	    {
	      type computer;
	      bmipsInt = 0 mips;
	      bmipsFP = 0 mips;
	      bpower = 0 watt;
	      dcost = 0 dollar;
	    }
	  
	} /* computers */
      
      Component Sensors
	{
	  Component CameraPanoSphere // FovHneed, FovVneed, freq are variables
	    {
	      type insulated;
	      resolution = 1024 * 1024 * 3 byte;
	      FovH = 360 degree;
	      FovV = 270 degree;
	      N = (int) FovHneed*FovVneed/(FovH*FovV);
	      freq = 1 Hz;
	      area = 20 cm * 20 cm * 5 * N;
	      mass = N * 3 kg;
	      /* Typecasting to (watt) does not work - FIx, hacked here ? */
	      power = (8 watt + 10 watt*(freq/14 Hz)) * N; /* OK? */
	      cost = 1.0 megadollar + 1.0 megadollar*N;
	      Component compression  : Fractal{}
	      // Syntax above OK?
	      mipsDint = N * resolution * freq * compression.computing;
	      downlink = N * resolution * freq / compression.ratio;
	      resolution <= 1024*1024*3 byte;
	      freq <= 14 Hz;
	    }
	  
	  Component CameraColor
	    {
	      resolution = 256 * 256 * 3 byte;
	      FovH = 90 degree;
	      FovV = 90 degree;
	      N = FovHneed*FovVneed/(FovH*FovV);
	      freq = 1 Hz;
	      area = 20 cm * 20 cm * 5 * N;
	      mass = N * 1 kg;
	      power =  N * (5 watt + 10 watt*(freq/30 Hz)); /* OK? */
	      cost = 1 + 1*N;
	      Component compression  : Fractal{}
	      mipsDint = N * resolution * freq * compression.computing;
	      downlink = N * resolution * freq / compression.ratio;
	      resolution <= 512*512*3 byte;
	      freq <= 30 Hz;
	    }
	}
      
      Component Software 
	{
	  Component Compression /* ratio for visually clean copy */
	    {
	      Component JPEG         
		{ratio =  10; computing = 445 Mips/(512*512*3 byte)/(30 Hz);}
	      Component MPEG2        
		{ratio =  70; computing = 445 Mips/(512*512*3 byte)/(30 Hz);}
	      Component Wavelet
		{
		  ratio = if(freq >= 5. hz,200,60); 
		  computing = if(freq>=5. hz,178 Mips/(512*512*3 byte)/(30 Hz),
				 178 Mips/(512*512*3 byte)/(30 Hz));
		}
	      Component Fractal
		{
		  ratio =  if(freq >=5. hz, 600,90); 
		  computing = if(freq >=5. hz, 178 Mips/(512*512*3 byte)/(30 Hz),
				 267 Mips/(512*512*3 byte)/(30 Hz));
		}
	      Component FractalScience
		{ratio =  40; computing = 267 Mips/(512*512*3 byte)/(30 Hz);}
	    }
	  Component Telemetry
	    {
	      Component MinimalTelemetry
		{
		  freq = 10 Hz; /* general rate people need info at */
		  packetsize = 4*20 byte; /* 20 doubles of state info */
		  downlink = freq * packetsize;
		}
	      Component MediumTelemetry
		{
		  freq = 10 Hz; /* general rate people need info at */
		  packetsize = 4*200 byte; 
		  /* 200 doubles of all sorts of state info */
		  downlink = freq * packetsize;
		}
	    } /* telemetry */
	  
	} /* software */
      
      Component SciencePackages
	{
	  Component SciencePackageLarge
	    {
	      Component Drill  : CryoDrillSATM{}
	      Component SampleAnalyser  : REGA{}
	      Component PanTilt  : LargePanTilt{}
	      Component IRCamera  : IRCameraClementine{}
	      Component IRLaser  : IRLaserClementine{}
	      Component HiResCamera  : HiResCameraClementine{}
	      Component Strobe  : StrobeGeneric{}
	      mipsGfp = 10 mips;
	    } /* SciencePackageLarge */
	  
	  Component SciencePackageSmall
	    {
	      Component HiResCamera  : HiResCameraClementine{}
	      Component PanTilt  : SmallPanTilt{}
	      Component Strobe  : StrobeGeneric{}
	    }
	  
	  Component SciencePackageNone
	    {
	      mass = 0.0 kg;
	      cost = 0.0 dollar;
	      power = 0.0 watt;
	    }
	  
	  Component CryoDrillSATM
	    {
	      type sensor, insulated;
	      cost = 2;
	      mass = 4 kg;
	      area = 1 m * 0.2 m * 4;
	      distrib power = distribution(max=15,avg=0,min=0) watt;
	      Component CryoSensors{type sensor; N = 8; /* 4 sensors, 4 actu. */ }
	      Component IRFiber {}
	      Component IRTransciever {}
	      Component Microscope {}
	    }
	  Component REGA
	    {
	      type sensor; /* already insulated */
	      cost = 0.5;
	      mass = 4 kg;
	      distrib power = distribution(max=20,avg=0,min=0) watt;
	    }
	  Component LargePanTilt /* mass est from Scheaffer mag numbers */
	    {
	      type insulated, sensor;
	      cost = 1;
	      mass = 3 kg; /* for two cameras, strobe and laser payload */
	      area = 20 cm * 30 cm * 5;
	      distrib power = distribution(max=10,avg=2,min=0) watt;
	      mipsGfp = 100 * 100 [byte/s];
	    }
	  Component SmallPanTilt /* mass est from Scheaffer mag numbers */
	    {
	      type insulated, sensor;
	      cost = 1;
	      mass = 1 kg; /* for two cameras, strobe and laser payload */
	      area = 20 cm * 15 cm * 5;
	      distrib power = distribution(max=6,avg=1,min=0) watt;
	      mipsGfp = 100 * 100 [byte/s];
	    }
	  Component IRCameraClementine
	    {
	      type insulated, sensor;
	      resolution = 256*256 byte;
	      freq = 0.5 Hz;
	      cost = 0.2;
	      mass = 2 kg;
	      area = 5 * 15 cm * 15 cm;
	      distrib power = distribution(max=5,avg=2,min=0) watt;
	      Component compression  : FractalScience{}
	      MipsDint = resolution*freq*compression.computing;
	      downlink = resolution*freq/compression.ratio;
	    }
	  Component IRLaserClementine
	    {
	      type insulated, sensor;
	      cost = 0.2;
	      mass = 2 kg;
	      area = 10 cm * 15 cm * 5;
	      distrib power = distribution(max=10,avg=2,min=0) watt;
	    }
	  Component HiResCameraClementine
	    {
	      type insulated, sensor;
	      resolution = 512*512 byte;
	      freq = 0.5 Hz;
	      cost = 0.5;
	      mass = 1.5 kg;
	      area = 15 cm * 15 cm * 5;
	      distrib power = distribution(max=8,avg=5,min=0) watt;
	      Component compression  : FractalScience{}
	      MipsDint = resolution*freq*compression.computing;
	      downlink = resolution*freq/compression.ratio;
	    }
	  Component StrobeGeneric
	    {
	      type insulated, sensor;
	      cost = 0.5;
	      mass = 1 kg;
	      area = 10 cm * 10 cm * 5;
	      distrib power = distribution(max=10,avg=4,min=0) watt;
	    }
	} /* SciencePackages */
    } /* Parts */
  
Component Missions
    {
      Component IceDiscovery
	/* Baseline  Lunar Ice Discovery Mission */
	{
	  duration = 30 day;
	  distrib speed = distribution(min=0.20, max=0.4, avg=0.25) [m/s];
	  distance = duration*get_min(Speed);
	  
	  latitude = 88 degree; // lat of mission
	  Tmin = planet.Tmin(latitude);  // Bug in degC to degK
	  Tmax = planet.Tmax(latitude);
	  
	  slopemax = 25 degree;
	  tiltNormal = 2 degree;
	  /* deviation of mission from landing g vector over mission lifetime
	     // how much does the normal tilt over the excursion geometry */
	  
	  // I suspect this metric is simplistic and inadequate, maybe harmful ?
	  obstacleMax = 0.15 m;
	  
	  // find a way to relate to battery over all missions ?
	  SortieDuration = 20 hour; /* hours */
	  NSorties = 10; /* max deep sorties */
	  energyNeed = SortieDuration * old Rover.power * 1.5; // needed in all missions
	  Night = 0.0 sec; // Needed in all missions 
	  
	  modify component Imager { FovHneed = 360 degree; FovVneed = 269 degree; freq = 1 Hz;}

	  modify component Receiver { D = 34 m;} // Earth Dish

	  modify component Rover { Component SciencePackage : SciencePackageLarge{} }
	  modify component Rover.Comm { Component transmitter : AnnularLink{}}
	  
	  modify component Rover.PowerSystem { Component PowerSource: solararray{}}
	}
      
      Component Edutainment
	{
	  distance = 2000 km;
	  duration = 2 year;
	  distrib speed = distribution(min=0.20, max = 0.5, avg = 0.4) [m/s];
	  
	  latitude = 0 degree;
	  Tmin =  planet.Tmin(latitude);
	  Tmax =  planet.Tmax(latitude); 
	  
	  slopemax = 25 degree;
	  tiltNormal = 30 degree;
	  obstacleMax = 0.15 m; // nothing, need a better obstacle metric
	  
	  NightPower = thermal.nightPowerNeed - sum(key NightPowerOut of rover);
	  energyNeed = NightPower * Night;
	  Night = planet.night; // How long do i have to live in the night
	  
	    modify component Imager {FovHneed = 360 degree; FovVneed = 269 degree; freq = 14 Hz;}

	  modify component Receiver { D = 34 m;} // Earth Dish

	  modify component Rover { Component SciencePackage : SciencePackageNone{}}
	  modify component Rover.Comm { Component transmitter : PhasedArray{}}

	  modify component Rover.PowerSystem { Component PowerSource: RTG{}}
	}
      
      Component ISRU /* pole vs equator */
	{
	  duration = 2 year;
	  distrib speed = distribution(min=0.20, max = 0.5, avg = 0.4) [m/s];
	  distance = duration*speed;
	  
	  latitude = 90;
	  Tmin = planet.Tmin(latitude);
	  Tmax = Planet.Tmax(latitude);
	  
	  slopemax = 15 degree;
	  tiltNormal = 1 degree;
	  ObstacleMax = 0.05 m;
	  
	  energyNeed = 0.0 joule;
	  Night = 0.0 sec; // no night as it hides in a warm barn
	 
	  modify component Imager {FovHneed = 360 degree; FovVneed = 269 degree; freq = 2 Hz;}
	  modify component Receiver { D = 22 m;} // Earth Dish
	  
	  modify component Rover { Component SciencePackage : SciencePackageNone{}}
	  modify component Rover.Comm { Component transmitter : PhasedArray{}}

	  modify component Rover.PowerSystem { Component PowerSource: SolarArray{}}
	}
      
      Component Millenium
	{
	  duration = 14 day;
	  distrib speed = distribution(min=0.10, max = 0.20, avg = 0.10);
	  distance = 5 km;
	  
	  latitude = 0 deg;
	  Tmin = Planet.Tdusk(latitude);
	  Tmax = Planet.Tmax(latitude);
	  
	  slopemax = 15 degree;
	  tiltNormal = 1 degree;
	  obstacleMax = 0.02 m;
	  
	  energyNeed = 0.0 joule;
	  Night = 0.0 sec;

	  modify component Imager {FovHneed = 90 degree; FovVneed = 89 degree; freq = 30 Hz;}
	  modify component Receiver { D = 34 m;} // Earth Dish
	  
	  modify component Rover { Component SciencePackage : SciencePackageNone{}}
	  modify component Rover.Comm { Component transmitter : PhasedArray{}}

	  modify component Rover.PowerSystem { Component PowerSource: solararray{}}
	}
      
    } /* Missions */
  
  Component Environment
    {
      component Moon
	{
	  RevPeriod = 28.45 day; /* in seconds */
	  Night = 14.5 day;
	  distrib EarthDist = Distribution(avg=384000, max=406697, min=356410) km; // actual avg may be 384000, plan for worst case
	  
	  cycle = 27.322 hour; /* approx max "day" hrs - sun cycle, see p57 LSB */
	  g = 9.81 [m/s2] / 6;
	  insolation = 1358 [watt/m2];
	  Tmax(latitude = 0 degree) = 130 degC - (185 degC)*latitude/(90 degree);
	  Tdusk(latitude = 0 degree) = 10 degC - (190 degC)*latitude/(90 degree);
	  Tmin(latitude = 0 degree) = -180 degC;
	}
      
      Component LunarSoil
	{ /* Lunar Soil LSS model, LSB, P526,8,.. */
	  /* Pressure sinkage characterstics under wheel load */
	  n = 1; /* regolith sinkage exponent */
	  Kc = 1382.32 [newton/m2]; /* cohesive modulus N/m^(n+1) */
	  Kphi = 816331.27 [newton/m3]; /* frictional modulus N/m^(n+2) */
	  
	  c = 172.79 [Newton/m2]; /* cohesion of regolith N/m^2 */
	  K = 0.01778 m; /* coeff of slip  meters */
	  phi = 35 degree; /* angle of internal friction/ repose 35 */
	  fr  = 0.1; /* coeffecient of rolling friction */
	  pressureMax = 6000 pascal; /* Pascals */
	  
	  /* needed by safeguarding & others */
	  mu_static = 0.3;
	  mu_rolling = 0.2; /* check these ? */
	  size = 0.005 m; /* is half mm reasonable soil irregularity ? */
	}
    }
  
  Component Analysis
    {
    } /* Analysis */
  
  Component SubSystems
    {
    } /* Subsystems */
  
} /* Options */



/***************************** Satellite Stuff *******************************/

Component OrbiterStage 
{
  dormant = true;
  Component Margin 
    {
      mass = 0.2*sum(key mass of OrbiterStage); /* same as power margin */
    }
  
  Component SRM 
    {
      Component Fuel{mass = 33.52 kg;}
      Component rocket 
	{
	  mass = 7.5 kg;
	  cost = 0.5;
	}
      Component interface 
	{
	  mass = 4 kg;
	  cost = 0.1;
	}
    } /* SRM */
  
  Component Satellite 
    {
      Component CylSolarArray
	{
	  mass = 2 kg;
	  cost = 0.1;
	}
      Component Gyro 
	{
	  mass = 8 kg;
	  power = 28;
	  cost = 0.5;
	}
      Component computerBoard 
	{
	  mass = 2 kg;
	  power = 10;
	  cost = 0.2;
	}
      Component cameraElectronics 
	{
	  mass = 3 kg;
	  distrib power = distribution(max=10,min=2);
	  cost = 0.4;
	}
      Component CommPlanet 
	{
	  mass = 3 kg;
	  power = 20;
	  cost = 0.5;
	}
      Component CommEarth 
	{
	  Component horn 
	    {
	      mass = 1 kg;
	    }
	  Component comm
	    {
	      mass = 5 kg;
	      power = 40;
	      cost = 1;
	    }
	  Component gimbal 
	    {
	      mass = 5 kg;
	      power = 5;
	      cost = 0.5;
	    }
	} /* CommEarth */
      Component battery 
	{
	  mass = 5.5 kg;
	  cost = 0.2;
	}
      Component PMAD 
	{
	  mass = 4 kg;
	  cost = 0.3;
	}
      
      Component structure 
	{
	  cost = 1;
	  Component support
	    {
	      mass = 0.2*sum(key mass of satellite);
	    }
	  Component bus
	    {
	      mass = 5 kg;
	    }
	  Component tankmount 
	    {
	      mass = 0.3*sum(key mass of coldgastank);
	    }
	} /* Structure */
      
      Component propulsion 
	{
	  cost = 0.5;
	  Component coldgastank 
	    {
	      Component coldgas {mass = 6 kg;}
	      Component tank{mass = 2 kg;}
	    } /* coldgastank */
	  Component nozzles {mass = 2 kg;}
	  Component actuators {mass = 2 kg;}
	} /* propulsion */
      
      Component mockup {cost = 0.6;}
    } /* satellite */
} /* orbiter stage */



Component CommSatellite
{ 
  dormant = true;
  /* Rover Comm, If there is satellite */
  mass = 3 kg;
  power = 10; /* for 1 W output */
  cost = 0.5 megadollar + backupEarth.cost;
  Component Transciever{}
  Component Modem{}
  Component Cable{}
  Component Antenna{}
  Component backupEarth
    {
      mass = 3 kg;
      distrib power = distribution(max=30,avg=0,min=0);
      cost = 1;
    }
} /* CommSatellite */

blah = wheel.rc( 0 radian, 1 newton );
blah2 = wheel.h( 0 radian, 1 newton );
blah3 = wheel.rb( 0 radian, 1 newton );
blaah = contactarea(wheelloadmax);
blaah2 = sinkage(wgt = wheelloadmax);


