Initial commit
[debian/openrocket] / src / net / sf / openrocket / aerodynamics / ExtendedISAModel.java
1 package net.sf.openrocket.aerodynamics;
2
3 import static net.sf.openrocket.aerodynamics.AtmosphericConditions.R;
4 import net.sf.openrocket.util.MathUtil;
5
6
7 /**
8  * An atmospheric temperature/pressure model based on the International Standard Atmosphere
9  * (ISA).  The no-argument constructor creates an {@link AtmosphericModel} that corresponds
10  * to the ISA model.  It is extended by the other constructors to allow defining a custom
11  * first layer.  The base temperature and pressure are as given, and all other values
12  * are calculated based on these.
13  * <p>
14  * TODO:  LOW:  Values at altitudes over 32km differ from standard results by ~5%.
15  * 
16  * @author Sampo Niskanen <sampo.niskanen@iki.fi>
17  */
18 public class ExtendedISAModel extends AtmosphericModel {
19         
20         public static final double STANDARD_TEMPERATURE = 288.15;
21         public static final double STANDARD_PRESSURE  = 101325;
22         
23         private static final double G = 9.80665;
24
25         private final double[] layer = {0, 11000, 20000, 32000, 47000, 51000, 71000, 84852};
26         private final double[] baseTemperature = {
27                         288.15, 216.65, 216.65, 228.65, 270.65, 270.65, 214.65, 186.95  
28         };
29         private final double[] basePressure = new double[layer.length]; 
30         
31
32         /**
33          * Construct the standard ISA model.
34          */
35         public ExtendedISAModel() {
36                 this(STANDARD_TEMPERATURE, STANDARD_PRESSURE);
37         }
38         
39         /**
40          * Construct an extended model with the given temperature and pressure at MSL.
41          * 
42          * @param temperature   the temperature at MSL.
43          * @param pressure              the pressure at MSL.
44          */
45         public ExtendedISAModel(double temperature, double pressure) {
46                 this(0, temperature, pressure);
47         }
48         
49         
50         /**
51          * Construct an extended model with the given temperature and pressure at the
52          * specified altitude.  Conditions below the given altitude cannot be calculated,
53          * and the values at the specified altitude will be returned instead.  The altitude
54          * must be lower than the altitude of the next ISA standard layer (below 11km).
55          * 
56          * @param altitude              the altitude of the measurements.
57          * @param temperature   the temperature.
58          * @param pressure              the pressure.
59          * @throws IllegalArgumentException  if the altitude exceeds the second layer boundary
60          *                                                                       of the ISA model (over 11km).
61          */
62         public ExtendedISAModel(double altitude, double temperature, double pressure) {
63                 if (altitude >= layer[1]) {
64                         throw new IllegalArgumentException("Too high first altitude: "+altitude);
65                 }
66                 
67                 layer[0] = altitude;
68                 baseTemperature[0] = temperature;
69                 basePressure[0] = pressure;
70                 
71                 for (int i=1; i < basePressure.length; i++) {
72                         basePressure[i] = getExactConditions(layer[i]-1).pressure;
73                 }
74         }
75         
76         
77         @Override
78         public AtmosphericConditions getExactConditions(double altitude) {
79                 altitude = MathUtil.clamp(altitude, layer[0], layer[layer.length-1]);
80                 int n;
81                 for (n=0; n < layer.length-1; n++) {
82                         if (layer[n+1] > altitude)
83                                 break;
84                 }
85
86                 double rate = (baseTemperature[n+1] - baseTemperature[n]) / (layer[n+1] - layer[n]);
87                 
88                 double t = baseTemperature[n] + (altitude - layer[n]) * rate;
89                 double p;
90                 if (Math.abs(rate) > 0.001) {
91                         p = basePressure[n] *
92                                 Math.pow(1 + (altitude-layer[n])*rate/baseTemperature[n], -G/(rate*R));
93                 } else {
94                         p = basePressure[n] *
95                                 Math.exp(-(altitude-layer[n])*G/(R*baseTemperature[n]));
96                 }
97                 
98                 return new AtmosphericConditions(t,p);
99         }
100
101         @Override
102         public double getMaxAltitude() {
103                 return layer[layer.length-1];
104         }
105         
106         
107         public static void main(String foo[]) {
108                 ExtendedISAModel model1 = new ExtendedISAModel();
109                 ExtendedISAModel model2 = new ExtendedISAModel(278.15,100000);
110                 
111                 for (double alt=0; alt < 80000; alt += 500) {
112                         AtmosphericConditions cond1 = model1.getConditions(alt);
113                         AtmosphericConditions cond2 = model2.getConditions(alt);
114                         
115                         AtmosphericConditions diff = new AtmosphericConditions();
116                         diff.pressure = (cond2.pressure - cond1.pressure)/cond1.pressure*100;
117                         diff.temperature = (cond2.temperature - cond1.temperature)/cond1.temperature*100;
118                         System.out.println("alt=" + alt + 
119                                         ": std:" + cond1 + " mod:" + cond2 + " diff:" + diff); 
120                 }
121         }
122
123 }