From abcb133abbc26703ec87dbc7c9799c432e276e29 Mon Sep 17 00:00:00 2001
From: Holger Niemann <holger.niemann@ipp.mpg.de>
Date: Wed, 22 Nov 2017 09:13:34 +0100
Subject: [PATCH] implement png download, multithread download, update FOV

---
 downloadversionIRdata.py | 205 ++++++++++++++++++++++++++++-----------
 1 file changed, 150 insertions(+), 55 deletions(-)

diff --git a/downloadversionIRdata.py b/downloadversionIRdata.py
index aec3ee9..d92d3ff 100644
--- a/downloadversionIRdata.py
+++ b/downloadversionIRdata.py
@@ -10,6 +10,10 @@ import W7Xrest.read_restdb as AKF_1
 import datetime
 import urllib
 import json
+import holn.archivedb as AKF_2
+from PIL import Image
+from io import BytesIO
+import multiprocessing
 
 archivepath="http://archive-webapi.ipp-hgw.mpg.de/Test/raw/W7X/"
 
@@ -437,6 +441,115 @@ def download_raw_images_by_times(port,starttime,stoptime,version=1,intervalSize=
         except urllib.error.URLError as e:
             print(e)
             return False, 0,-1
+
+
+def download_raw_images_by_program_via_png(port,program,time_s=0,version=1,threads=1):
+    prog=AKF_1.get_program_from_PID(program)
+    if prog[0]:
+        starttime=prog[1]['trigger']['1'][0]
+        stoptime=prog[1]['trigger']['6'][0]
+        success=True
+        OP=get_OP_by_time(starttime)
+        Cam=portcamdict[OP]['AEF'+str(port)]
+        if Cam.split("_")[0]=="Infratec":#camera=="INFRATEC" or camera=="infratec" or camera=="Infratec":
+            larchivepath="Test/raw/W7X/"+"QRT_INFRATEC/"+"AEF"+str(port)+"_raw_DATASTREAM/V"+str(version)+"/0/raw"
+        elif Cam.split("_")[0]=="IRCam":#camera=="IRCAM" or camera=="IRcam" or camera=="ircam":
+            larchivepath="Test/raw/W7X/"+"QRT_IRCAM/"+"AEF"+str(port)+"_raw_DATASTREAM/V"+str(version)+"/0/raw"
+        else:
+            raise Exception("Port number does not fit the known cameras")
+        stdate=datetime.datetime.utcfromtimestamp((starttime-100)/1e9)
+        stdate=stdate.isoformat()
+        if time_s==0:
+            enddate=datetime.datetime.utcfromtimestamp(stoptime/1e9)        
+            enddate=enddate.isoformat()
+        else:
+            enddate=datetime.datetime.utcfromtimestamp((starttime)/1e9+time_s)  
+            enddate=enddate.isoformat()
+            #"2017-11-15 08:00:00"
+        times=AKF_2.get_all_time_intervals(larchivepath,stdate.replace("T"," "),enddate.replace("T"," "))
+        time=[]
+        images=[]
+        lnt=len(times)
+        if threads==1:
+            for i in range(lnt):
+                ele=times[lnt-1-i]
+                imag=download_last_raw_image_by_time(port,ele[0]-10,ele[0]+10)
+                if imag[0]:
+                    time.append(ele[0])
+                    images.append(imag[1])
+                else:
+                    success=False
+            return success,np.array(time),np.array(images,dtype=np.uint16)
+        else:
+            tim=[]
+            for i in range(lnt):
+                tim.append(times[lnt-1-i][0])
+            intervalls=[]
+            intervalSize=int(lnt/threads)
+            for i in range(threads):
+                intervalls.append(int(i*intervalSize))
+            intervalls.append(lnt)
+            jobs = []
+            out_q=multiprocessing.Queue()
+            for i in range(threads):
+                print("Start Thread ",i+1)                
+                p = multiprocessing.Process(target=download_raw_images_png_by_times_thread, args=(port,tim[intervalls[i]:intervalls[i+1]],out_q,i,version,))
+                jobs.append(p)
+                p.start()
+            resultdict = []
+            for i in range(threads):
+                resultdict.append(out_q.get())
+            for p in jobs:
+                p.join()
+            order=[]
+            for ele in resultdict:
+                order.append(ele[0])
+                if len(np.where(np.asarray(ele[1])==False)[0])>0:
+                    success=False
+            
+            times=np.array(resultdict[order.index(0)][2])
+            images=np.array(resultdict[order.index(0)][3])
+            for i in range(threads-1):
+                images=np.append(images,resultdict[order.index(i+1)][3],axis=0)
+                times=np.append(times,resultdict[order.index(i+1)][2])
+            del resultdict
+            return success,times,images
+        
+
+def download_raw_images_png_by_times_thread(port,times,out_q,threadnumber,version=1):
+    images=[]
+    time=[]
+    successes=[]
+    for i in times:
+        imag=download_last_raw_image_by_time(port,i-10,i+10)
+        if imag[0]:
+            images.append(imag[1])
+            time.append(i)
+            successes.append(True)
+        else:
+            successes.append(False)
+    out_q.put([threadnumber,successes,time,images])        
+#    return success,np.array(time),np.array(images,dtype=np.uint16)
+
+def download_last_raw_image_by_time(port,starttime,stoptime,version=1):
+    OP=get_OP_by_time(starttime)
+    Cam=portcamdict[OP]['AEF'+str(port)]
+    if Cam.split("_")[0]=="Infratec":#camera=="INFRATEC" or camera=="infratec" or camera=="Infratec":
+        larchivepath=archivepath+"QRT_INFRATEC/"+"AEF"+str(port)+"_raw_DATASTREAM/V"+str(version)+"/0/raw"
+    elif Cam.split("_")[0]=="IRCam":#camera=="IRCAM" or camera=="IRcam" or camera=="ircam":
+        larchivepath=archivepath+"QRT_IRCAM/"+"AEF"+str(port)+"_raw_DATASTREAM/V"+str(version)+"/0/raw"
+    else:
+        raise Exception("camera unknown, stopping here")
+    try:
+        res = urllib.request.urlopen(larchivepath+"/_signal.png?from="+str(starttime-10)+"&upto="+str(stoptime))
+        img = Image.open(BytesIO(res.read()))
+        res.close()
+#        pixelarray = np.array(img.getdata()).reshape(img.size[1], img.size[0])
+        pixelarray = np.array(img,dtype=np.uint16)#.swapaxes(0,1)
+        return True, pixelarray
+    except urllib.error.URLError as e:
+        print(e)
+        return False, -1
     
 def download_raw_parlog_by_program(port,program,version=1):
     prog=AKF_1.get_program_from_PID(program)
@@ -990,10 +1103,8 @@ def get_temp_from_raw_by_program_V2(portnr,program,time_s=0,emi=0.8,tubepart=0,v
                             raise Exception
                         else:
                             gain=NUC_DL[1][0]
-                            offset=NUC_DL[1][1]
-                            
-                            print(datetime.datetime.now(),"Start download of raw images")
-                            
+                            offset=NUC_DL[1][1]                            
+                            print(datetime.datetime.now(),"Start download of raw images")                            
                             if (stoptime-starttime)/intervalSize>1:        
                                 nrinterv=int(np.ceil((stoptime-starttime)/intervalSize))
                                 print("timewindow to large, putting it into smaller fractions ("+str(nrinterv)+")")
@@ -1196,7 +1307,7 @@ def get_calib_data(port,program,emissivity=0.8,Temp_V=2,version=1):
                             offset=NUC_DL[1][1]
                             badpixels=NUC_DL[1][3]
                             if np.max(badpixels)==0:
-                                badpixels=find_badpixels(gain,offset)
+                                badpixels=find_badpixels(port,gain,offset)
                     else:
                         gain=0
                         offset=0
@@ -1206,7 +1317,7 @@ def get_calib_data(port,program,emissivity=0.8,Temp_V=2,version=1):
                 raise Exception("no LUT found")
         else:
             raise Exception("no exposure time found")
-        return background,LUT,refT,gain,offset
+        return background,LUT,refT,gain,offset,badpixels
     else:
         print("Program was not found")
         return 0,0,0,0,0
@@ -1324,48 +1435,22 @@ def restore_pixels(frame, bad_pixel):
             print('ERROR: no adjacent pixel found!')
     return resframe
 
-#def get_average_background_circle(port,image):
-#    try:
-#        points=valid_background_circle[port]
-#        dummy=[]
-#        for y in range(len(image)):
-#            """
-#            Sekante erreichnen für beide Kreise, Gerade ax+by=c, sonderfall: y=c, a=0,b=1
-#            Kreis definiert durch (x-x0)²+(y-y0)²=r²
-#            loesung:  d=c-a*x0-b*y0=c-y0
-#                x1,2=x0+(ad+-b*sqrt(r²(a²+b²)-d²))/(a²+b²)=x0+-sqrt(r²-d²)
-#                y1,2=y0+(bd-+a*sqrt(r²(a²+b²)-d²))/(a²+b²)=y0+d=c=y
-#            """
-#            y0=points[1]
-#            x0=points[0]
-#            r1=points[2]
-#            r2=points[3]
-#            xs=[x0+np.sqrt(r1**2-(y-y0)**2),x0+np.sqrt(r2**2-(y-y0)**2),x0-np.sqrt(r1**2-(y-y0)**2),x0-np.sqrt(r2**2-(y-y0)**2)]
-#            xs.sort()
-#            for i in range(4):
-#                if xs[i]<0:
-#                    xs[i]=0
-#                elif xs[i]>np.shape(image)[1]:
-#                    xs[i]=np.shape(image)[1]
-#            for xi in range(int(xs[0]+1),int(xs[1])):
-#                dummy.append(image[y][xi])
-##                back[y][xi]=0
-#            for xi in range(int(xs[2]+1),int(xs[3])):
-#                dummy.append(image[y][xi])
-##                back[y][xi]=0
-#        return np.average(dummy)
-#    except Exception as E:
-#        print(E)
-#        return 0
 
 def make_FOV_mask(port):
     points=valid_FOV_circle[port]
+    """
+    Sekante erreichnen für Kreis, Gerade ax+by=c, sonderfall: y=c, a=0,b=1
+    Kreis definiert durch (x-x0)²+(y-y0)²=r²
+    loesung:  d=c-a*x0-b*y0=c-y0
+        x1,2=x0+(ad+-b*sqrt(r²(a²+b²)-d²))/(a²+b²)=x0+-sqrt(r²-d²)
+        y1,2=y0+(bd-+a*sqrt(r²(a²+b²)-d²))/(a²+b²)=y0+d=c=y
+    """
     y0=points[1]
     x0=points[0]
     r1=points[2]
-    da,time,back=download_background_by_program(port,"20171115.034",9)
+    da,time,back=download_background_by_program(port,"20171109.045",9)
     fig = plt.figure()
-    plt.imshow(back,vmin=np.average(back)-500,vmax=np.average(back)+1500)
+    plt.imshow(back,vmin=np.average(back)-200,vmax=np.average(back)+500)
     inner_circle = mlt.patches.Circle((x0,y0), r1,color = 'r', fill = False)
     ax = fig.gca()
     ax.add_artist(inner_circle)
@@ -1401,7 +1486,15 @@ def get_average_background_recangle(port,image):
 
 valid_FOV_circle = {
     10: [562, 364, 550],#, 600],
-    11: [502, 404, 520]#, 570]
+    11: [502, 404, 520],#, 570]
+    20: [535, 374, 510],
+    21: [548, 404, 505],
+    30: [517, 389, 520],
+    31: [562, 364, 520],
+    40: [542, 394, 520],
+    41: [522, 394, 520],
+    50: [640, 454, 620],
+    51: [512, 379, 510]
     }
 
 valid_background_rectangle = {
@@ -1414,25 +1507,27 @@ if __name__=='__main__':
     import matplotlib.pyplot as plt
     import matplotlib as mlt
     
-    make_FOV_mask(10)
+#    make_FOV_mask(31)
+    
 #    bla=get_average_background_recangle(11,back)
 #    print(bla)
-#    exist,time,images=download_raw_images_by_program(10,"20171115.034",0.1)
-#    background1,LUT,refT1,gain,offset=get_calib_data(10,"20171115.034",0.8,1)
+    background1,LUT,refT1,gain,offset,badpixels=get_calib_data(10,"20171115.034",0.8,1)
+    plt.figure()
+    plt.imshow(badpixels)
+    plt.figure()
+    plt.imshow(gain,vmin=0,vmax=10)
+    badpixels=np.zeros(np.shape(gain))
+    certain_bads=np.zeros(np.shape(gain))
+    gainmax=100
+    certain_bads+=(gain>gainmax)
+    certain_bads+=(gain<0)
+    possi_goodgain=gain*(1-certain_bads)
+    plt.hist(possi_goodgain.flatten(),100,(0,10))
 #    background2,LUT,refT2,gain,offset=get_calib_data(10,"20171115.034",0.8,2)
 #    temp1=apply_calib_on_raw(images,background1,LUT,refT1,gain,offset,True)
 #    temp2=apply_calib_on_raw(images,background2,LUT,refT1,gain,offset,False)
 #    val1,time1,temp1,trangeval1=get_temp_from_raw_by_program_V1(10,"20171115.034",0.9)
 #    val2,time2,temp2,trangeval2=get_temp_from_raw_by_program_V2(10,"20171115.034",0.9,1)
 #    val3,time3,temp3,trangeval3=get_temp_from_raw_by_program_V2(10,"20171115.034",0.9,0)
-#    plt.figure()
-#    plt.imshow(np.asarray(temp1[-1])-273.15,vmin=50,vmax=400)
-#    plt.colorbar()
-#    plt.figure()
-#    plt.imshow(np.asarray(temp2[-1])-273.15,vmin=50,vmax=400)
-#    plt.colorbar()
-#    plt.figure()
-#    plt.imshow(np.asarray(temp3[-1])-273.15,vmin=50,vmax=400)
-#    plt.colorbar()
-#    plt.show()
+    plt.show()
     
\ No newline at end of file
-- 
GitLab