# This scripts is giving a quick n dirty overview of single counter in given time span
# Patrick Kuehl, 17.06.2015
import time
#from pylab import *
from datetime import datetime
import matplotlib.pyplot as plt
import numpy as np
from matplotlib import dates
import time
import datetime as date
from datetime import datetime
from datetime import timedelta
import os
import warnings


def load_soho_ephin_level1(start_doy,end_doy,year,filename):
  startyear,endyear=year,year
  mean=4  #
  init=1
  for thisyear in range(startyear,endyear+1): 
    sdoy=start_doy
    edoy=end_doy
    for thisdoy in range(sdoy,edoy):
      print "Loading DoY %s of %s..." %(thisdoy,thisyear)
      if thisyear <2000: 
        thisyear2d=thisyear-1900
        prefix='eph'
      else: 
        thisyear2d=thisyear-2000
        prefix='epi'
      try:
        ephindata=np.loadtxt("/data/missions/soho/costep/level1/sci/%s/%s%02d%03d.sci" %(thisyear,prefix, thisyear2d,thisdoy))
        yearephin=ephindata[:,0]
        dayephin=ephindata[:,1]
        millisecephin=ephindata[:,2]
        g=ephindata[:,5]
        a=ephindata[:,6]+ephindata[:,7]+ephindata[:,8]+ephindata[:,9]+ephindata[:,10]+ephindata[:,11]
        b=ephindata[:,12]+ephindata[:,13]+ephindata[:,14]+ephindata[:,15]+ephindata[:,16]+ephindata[:,17]
        c=ephindata[:,18]
        d=ephindata[:,19]
        e=ephindata[:,20]
        f=ephindata[:,21]
        status=ephindata[:,-1]
        ringoff=np.zeros(len(status))
        for q in range(len(status)):
          binaries='{0:08b}'.format(int(status[q]))
          if int(binaries[-2]): ringoff[q]=1

        thisoutput=np.array(([yearephin,dayephin,millisecephin,g,a,b,c,d,e,f,ringoff])).T
        if init==0:
          output=np.vstack((output,thisoutput))
        if init==1:
          output=thisoutput
          init=0
      except: dummy=1
  #np.savetxt("data/%s" %filename, output, fmt='%i')
  return output

my_year=int(raw_input("Insert year of interest: ")) #2007
start_doy=int(raw_input("Insert start DoY: "))
end_doy=int(raw_input("Insert end DoY: "))

def fill_ringoff(tax):
  tax.fill_between(dateephin, tax.get_ylim()[0], tax.get_ylim()[1], where=ringoff>0, facecolor='red', interpolate=False,alpha=0.3)


for year in [my_year]: #1995
  output=load_soho_ephin_level1(start_doy,end_doy,year,"%s.dat" %year)
  yearephin,dayephin,millisecephin,g,a,b,c,d,e,f,ringoff=output[:,0],output[:,1],output[:,2], output[:,3],output[:,4],output[:,5],output[:,6],output[:,7],output[:,8],output[:,9],output[:,10]
  dateephin=[(datetime(int(yearephin[date]),1,1,0)+timedelta(int(dayephin[date])-1,0,0,int(millisecephin[date]))) for date in range(len(yearephin))]
  my_fmt="k-"

  plt.figure()
  p1=plt.subplot(711)
  plt.setp( p1.get_xticklabels(), visible=False)
  plt.plot_date(dateephin,g,my_fmt,label='G')
  plt.yscale('log')
  plt.ylabel("G")
  plt.xlim([datetime(year,1,1),datetime(year,12,31)])
  plt.ylim([10**3,10**7])
  p1.set_xticklabels([])
  plt.title("single count rates")
  fill_ringoff(plt.gca())

  p2=plt.subplot(712, sharex=p1)
  plt.setp( p2.get_xticklabels(), visible=False)
  plt.plot_date(dateephin,a,my_fmt,label='A')
  plt.yscale('log')
  plt.ylabel("A")
  plt.xlim([datetime(year,1,1),datetime(year,12,31)])
  plt.ylim([10**3,10**7])
  p2.set_xticklabels([])
  fill_ringoff(plt.gca())

  p3=plt.subplot(713, sharex=p1)
  plt.setp( p3.get_xticklabels(), visible=False)
  plt.plot_date(dateephin,b,my_fmt,label='B')
  plt.yscale('log')
  plt.ylabel("B")
  plt.xlim([datetime(year,1,1),datetime(year,12,31)])
  plt.ylim([10**3,10**7])
  p3.set_xticklabels([])
  fill_ringoff(plt.gca())

  p4=plt.subplot(714, sharex=p1)
  plt.setp( p4.get_xticklabels(), visible=False)
  plt.plot_date(dateephin,c,my_fmt,label='C')
  plt.yscale('log')
  plt.ylabel("C")
  plt.xlim([datetime(year,1,1),datetime(year,12,31)])
  plt.ylim([10**3,10**7])
  p4.set_xticklabels([])
  fill_ringoff(plt.gca())

  p5=plt.subplot(715, sharex=p1)
  plt.setp( p5.get_xticklabels(), visible=False)
  plt.plot_date(dateephin,d,my_fmt,label='D')
  plt.yscale('log')
  plt.ylabel("D")
  plt.xlim([datetime(year,1,1),datetime(year,12,31)])
  plt.ylim([10**3,10**7])
  p5.set_xticklabels([])
  fill_ringoff(plt.gca())

  p6=plt.subplot(716, sharex=p1)
  plt.setp( p6.get_xticklabels(), visible=False)
  plt.plot_date(dateephin,e,my_fmt,label='E')
  plt.yscale('log')
  plt.ylabel("E")
  plt.xlim([datetime(year,1,1),datetime(year,12,31)])
  plt.ylim([10**3,10**7])
  p6.set_xticklabels([])
  fill_ringoff(plt.gca())

  p7=plt.subplot(717, sharex=p1)
  plt.plot_date(dateephin,f,my_fmt,label='F')
  plt.yscale('log')
  plt.ylabel("F")
  #plt.xlim([datetime(year,1,1),datetime(year,12,31)])
  plt.ylim([10**3,10**7])
  hfmt=dates.DateFormatter('%d.%b \n %Y \n %H:%M')  #'%b \n %Y')  '%H:%M \n %D'
  p7.xaxis.set_major_formatter(hfmt)
  fill_ringoff(plt.gca())
  plt.show()
  #plt.savefig("counts_%s.pdf" %year)
  #plt.savefig(pdf, format='pdf') # note the format='pdf' argument!