//Arduino IDE 1.0.5
//www.basicairdata.eu Barometric altitude VTI SCP1000 sensor, JLJ
//SPI handling by http://forum.arduino.cc/ site http://forum.arduino.cc/index.php?topic=19311.45;wap2
//Infinite reading loop
//Pressure data smoothed with a simple moving average filter, non weighted mean of current value of the last two samples
// define spi bus pins
#define SLAVESELECT 10  //CSB
#define SPICLOCK 13     //SCK
#define DATAOUT 11	//MOSI
#define DATAIN 12	 //MISO
#define UBLB(a,b)  ( ( (a) << 8) | (b) )
#define UBLB19(a,b) ( ( (a) << 16 ) | (b) )

//Addresses
#define REVID 0x00	//ASIC Revision Number
#define OPSTATUS 0x04   //Operation Status
#define STATUS 0x07     //ASIC Status
#define START 0x0A      //Constant Readings
#define PRESSURE 0x1F   //Pressure 3 MSB
#define PRESSURE_LSB 0x20 //Pressure 16 LSB
#define TEMP 0x21       //16 bit temp
#define FIELD_ALTITUDE 53 //calibration height

char rev_in_byte;  //Read value from sensor register  
int temp_in;//Read temp value from sensor register
long base_pascal=99700;  //Pressure at sea level, from outside source
unsigned long pressure_lsb;  //Least significant byte
unsigned long pressure_msb; //Most significant byte
unsigned long pressure=0; //Pressure reading from sensor
unsigned long pressurem_one; //One iteration old pressure reading value 
unsigned long pressurem_two;//Two iterations old pressure reading value
unsigned long filtered; //Smoother output
unsigned long init_pressure; //Calculated initial sea level pressure from calibration
float init_temp;//Calculated initial sea level temperature from calibration
boolean init_loop=true; //Set the initial context
float RE=6137;//Local Earth radius

void setup()
{
  byte clr;
  pinMode(DATAOUT, OUTPUT);
  pinMode(DATAIN, INPUT);
  pinMode(SPICLOCK,OUTPUT);
  pinMode(SLAVESELECT,OUTPUT);
  digitalWrite(SLAVESELECT,HIGH); //disable device  

  SPCR = B01010011; //MPIE=0, SPE=1 (on), DORD=0 (MSB first), MSTR=1 (master), CPOL=0 (clock idle when low), CPHA=0 (samples MOSI on rising edge), SPR1=0 & SPR0=0 (500kHz)
  clr=SPSR;
  clr=SPDR;
  delay(10);
  Serial.begin(9600);
  delay(600);

  Serial.println("Continuous reading starts");
  write_register(0x03,0x09); //High speed acquisition mode start, max 9Hz, 15 bits
}

void loop()
{

  //Initialitation loop, it is executed only one time or when init_loop is set
  if (init_loop==true) {
    delay(300);
    for (int i=0;i<9;i++)
    {
      pressurem_two=pressurem_one;
      pressurem_one=pressure;
      rev_in_byte = read_register(REVID);
      pressure_msb = read_register(PRESSURE);
      pressure_msb &= B00000111;
      pressure_lsb = read_register16(PRESSURE_LSB);
      pressure = UBLB19(pressure_msb, pressure_lsb);
      pressure /= 4;
      init_pressure=pressure/3+pressurem_one/3+pressurem_two/3;
    }


    init_temp = read_register16(TEMP);
    init_temp = init_temp / 20;
    init_temp=init_temp +0.0065*(RE*FIELD_ALTITUDE/(RE+FIELD_ALTITUDE));
    Serial.print("Initialization sea level Temperature C ");
    Serial.println(init_temp ,1);
    init_pressure=init_pressure/(pow((1-0.0065*FIELD_ALTITUDE/(init_temp+273.15)),(5.2557862)));
    Serial.print("Initialization sea level pressure  Pa ");
    Serial.println(init_pressure, DEC);
    Serial.print("Sea level pressure from weather service Pa ");
    Serial.println(base_pascal, DEC);
    init_loop=false;  
  }

  //Moving average smoother, deep 2 past samples
  pressurem_two=pressurem_one;
  pressurem_one=pressure;
  rev_in_byte = read_register(REVID);
  pressure_msb = read_register(PRESSURE);
  pressure_msb &= B00000111;
  pressure_lsb = read_register16(PRESSURE_LSB);
  pressure = UBLB19(pressure_msb, pressure_lsb);
  pressure /= 4;
  filtered=pressure/3+pressurem_one/3+pressurem_two/3;  
  Serial.print("Raw pressure Pa ");
  Serial.println(pressure, DEC);
  Serial.print("Filtered pressure Pa ");
  Serial.println(filtered, DEC);
  temp_in = read_register16(TEMP);
  temp_in = temp_in / 20;
  Serial.print("Temperature C ");
  Serial.println(temp_in , DEC);
  Serial.print("Barometric altitude according, ISA m ");
  Serial.println(standardday(filtered),DEC);
  Serial.print("Barometric altitude corrected for pressure at sea level, m ");
  Serial.println(correctedstandardpressure(filtered),DEC);
  Serial.print("Barometric altitude corrected, current altitude and temperature known m ");
  Serial.println(correctedaltitude(filtered),DEC);
  //main delay
  delay(300);
}

char spi_transfer(volatile char data)
{
  SPDR = data;			  // Start the transmission
  while (!(SPSR & (1<<SPIF)))     // Wait for the end of the transmission
  {
  };
  return SPDR;			  // return the received byte
}
char read_register(char register_name)
{
  char in_byte;
  register_name <<= 2;
  register_name &= B11111100; //Read command

  digitalWrite(SLAVESELECT,LOW); //Select SPI Device
  spi_transfer(register_name); //Write byte to device
  in_byte = spi_transfer(0x00); //Send nothing, but we should get back the register value
  digitalWrite(SLAVESELECT,HIGH);
  delay(10);
  return(in_byte);

}
float read_register16(char register_name)
{
  byte in_byte1;
  byte in_byte2;
  float in_word;

  register_name <<= 2;
  register_name &= B11111100; //Read command

  digitalWrite(SLAVESELECT,LOW); //Select SPI Device
  spi_transfer(register_name); //Write byte to device
  in_byte1 = spi_transfer(0x00);    
  in_byte2 = spi_transfer(0x00);
  digitalWrite(SLAVESELECT,HIGH);
  in_word = UBLB(in_byte1,in_byte2);
  return(in_word);
}

void write_register(char register_name, char register_value)
{
  register_name <<= 2;
  register_name |= B00000010; //Write command

  digitalWrite(SLAVESELECT,LOW); //Select SPI device
  spi_transfer(register_name); //Send register location
  spi_transfer(register_value); //Send value to record into register
  digitalWrite(SLAVESELECT,HIGH);
}
//Function return altitude in standard day condition
float standardday(unsigned long rawpressure )
{
  float rawpressure_psi; //Pressure reading in psi, to follow references
  float base_psi;
  float h; //Geopotential altitude
  float H; //Geometric altitude
  rawpressure_psi=rawpressure/133.3223684211/25.4; //
  base_psi=101325/133.3223684211/25.4;
  h=0.3048*(pow(base_psi,0.190255)-pow(rawpressure_psi,0.190255))/0.000013125214; 
  H=h*RE/(RE-h);
  return(H);
}
//Function returns altitude corrected for sea level pressure deviation
float correctedstandardpressure(unsigned long rawpressure )
{
  float rawpressure_psi; //Pressure reading in psi, to follow references
  float base_psi;
  float h; //Geopotential altitude
  float H; //Geometric altitude
  rawpressure_psi=rawpressure/133.3223684211/25.4; //
  base_psi=base_pascal/133.3223684211/25.4;
  h=0.3048*(pow(base_psi,0.190255)-pow(rawpressure_psi,0.190255))/0.000013125214;
  H=h*RE/(RE-h);
  return(H);
}

//Function returns altitude corrected for pressure and temperature deviations
float correctedaltitude(unsigned long rawpressure)
{
  float rawpressure_psi; //Pressure reading in psi, according to references
  float base_psi;
  float h; //Geopotential altitude
  float H; //Geometric altitude
  //init_temp, init_pressure  
  rawpressure_psi=rawpressure/133.3223684211/25.4; //
  base_psi=init_pressure/133.3223684211/25.4;
  //Imperial units
  // h=0.3048*(pow(base_psi,0.190255)-pow(rawpressure_psi,0.190255))/0.000013125214*(init_temp+273.15)/(15+273.15); //Hypsometric formula
  //SI units
  h=(pow((float(init_pressure)/float(rawpressure)),0.190255)-1)*(init_temp+273.15)/0.0065;
  H=h*RE/(RE-h);
  return(H);
}