search

leginon-org / leginon-redmine-archive

1 stars 0 forks source link

Allow a cs value of 0 when uploading images #2762

Open leginonbot opened 8 months ago

leginonbot commented 8 months ago

Author Name: Amber Herold (Amber Herold) Original Redmine Issue: 2762, https://emg.nysbc.org/redmine/issues/2762 Original Date: 2014-05-07 Original Assignee: Michael Cianfrocco

ref: http://emg.nysbc.org/redmine/boards/13/topics/2049?r=2053

The gui and python script should allow cs corrected images to be uploaded.

Merge to 3.0

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Amber Herold (Amber Herold) Original Date: 2014-05-07T16:03:36Z

OK Mike, Go ahead and give this a try. Let me know if you run into any issues.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Michael Cianfrocco (@mcianfrocco) Original Date: 2014-05-07T16:16:59Z

Hi Amber,

I was able to use the GUI to generate a command, but then I got a command line error: 

  File "/labdata/allab/michaelc/myami-3.0-trunk//appion/imageloader.py", line 581, in <module>
    imgLoop = ImageLoader()
  File "/labdata/allab/michaelc/myami-3.0-trunk//appion/imageloader.py", line 30, in __init__
    appionLoop2.AppionLoop.__init__(self)
  File "/labdata/allab/michaelc/myami-3.0-trunk/appion/appionlib/appionLoop2.py", line 31, in __init__
    appionScript.AppionScript.__init__(self)
  File "/labdata/allab/michaelc/myami-3.0-trunk/appion/appionlib/appionScript.py", line 84, in __init__
    self.checkConflicts()
  File "/labdata/allab/michaelc/myami-3.0-trunk//appion/imageloader.py", line 102, in checkConflicts
    apDisplay.printError("Cs value must be in mm (e.g., 2.0)")
  File "/labdata/allab/michaelc/myami-3.0-trunk/appion/appionlib/apDisplay.py", line 65, in printError
    raise Exception, colorString("\n *** FATAL ERROR ***\n"+text+"\n\a","red")
Exception: 
 *** FATAL ERROR ***
Cs value must be in mm (e.g., 2.0)

Exception AttributeError: "'NoneType' object has no attribute 'startswith'" in <bound method ImageLoader.__del__ of <__main__.ImageLoader object at 0x27e1cd0>> ignored
leginonbot commented 8 months ago

Original Redmine Comment Author Name: Amber Herold (Amber Herold) Original Date: 2014-05-07T17:04:03Z

I'm not getting this error on my end. Can you double check that you are at r18253?

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Amber Herold (Amber Herold) Original Date: 2014-05-08T15:23:00Z

Mike, lets continue this just in the issue at this point. Your last post:

Hi Amber,

I got around that problem I posted to the issues page, but now i'm getting an error:

 ==== Committing data to database ==== 
 ... Committing ctf parameters for 14may07z_14apr11z_A1_49k_1 to database
!!! WARNING: |def1| > |def2|, flipping defocus axes
 ... Reading image...
[CTF param] def1: 6.37e-07 | def2: 6.89e-07 | angle: -74.1 | ampcontr 0.15 | defratio 1.082
Traceback (most recent call last):
  File "/labdata/allab/michaelc/myami-3.0-trunk//appion/bin/ctfestimate.py", line 421, in <module>
    imgLoop.run()
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/appionLoop2.py", line 97, in run
    self.loopCommitToDatabase(imgdata)
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/appionLoop2.py", line 146, in loopCommitToDatabase
    return self.commitToDatabase(imgdata)
  File "/labdata/allab/michaelc/myami-3.0-trunk//appion/bin/ctfestimate.py", line 329, in commitToDatabase
    ctfinsert.validateAndInsertCTFData(imgdata, self.ctfvalues, self.ctfrun, self.params['rundir'])
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctfinsert.py", line 45, in validateAndInsertCTFData
    ctfvalues = runCTFdisplayTools(imgdata, ctfvalues, opimagedir, fftpath, fftfreq)
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctfinsert.py", line 80, in runCTFdisplayTools
    ctfdisplaydict = ctfdisplay.makeCtfImages(imgdata, ctfvalues, fftpath, fftfreq)
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctfdisplay.py", line 1120, in makeCtfImages
    ctfdisplaydict = a.CTFpowerspec(imgdata, ctfdata, fftpath, fftfreq, twod=twod)
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctfdisplay.py", line 972, in CTFpowerspec
    self.convertDefociToConvention(ctfdata)
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctfdisplay.py", line 916, in convertDefociToConvention
    self.initfreq, self.cs, self.volts, self.ampcontrast)
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctftools.py", line 162, in defocusRatioToEllipseRatio
    cs, volts, ampcontrast, numzeros=1, zerotype="valleys")
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctftools.py", line 71, in getCtfExtrema
    radsq1 = (-b + root)/(2*a)
ZeroDivisionError: float division
It looks like Cs=0 results in an equation that tries to divide by zero.

Sorry to keep troubling you with this.

mike

Yep, I was a bit worried that might happen considering we were enforcing a min value of 0.1. I'm out of the office until Monday, and can see if there is a way around this issue at that time.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Michael Cianfrocco (@mcianfrocco) Original Date: 2014-05-08T15:36:56Z

It looks like Ace2 outputs a separate error, just to document it: 

 ... Angle astigmatism: 2.65 degrees
 ... Amplitude contrast: 50.00 percent
Final confidence: 0.533
 ==== Committing data to database ==== 
 ... Looking up session, 14may07z
 ... Committing ctf parameters for 14may07z_14apr11z_A1_49k_10 to database
!!! WARNING: atypical defocus #1 value -0.0000 microns (underfocus is positve)
!!! WARNING: Bad CTF values, insert but not create images
Traceback (most recent call last):
  File "/labdata/allab/michaelc/myami-3.0-trunk//appion/bin/pyace2.py", line 361, in <module>
    imgLoop.run()
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/appionLoop2.py", line 97, in run
    self.loopCommitToDatabase(imgdata)
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/appionLoop2.py", line 146, in loopCommitToDatabase
    return self.commitToDatabase(imgdata)
  File "/labdata/allab/michaelc/myami-3.0-trunk//appion/bin/pyace2.py", line 286, in commitToDatabase
    ctfinsert.validateAndInsertCTFData(imgdata, self.ctfvalues, self.ctfrundata, self.params['rundir'])
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctfinsert.py", line 71, in validateAndInsertCTFData
    ctfdb.printCtfData(ctfq)
  File "/labdata/allab/michaelc/myami-3.0-trunk/myami/lib/python2.6/site-packages/appionlib/apCtf/ctfdb.py", line 48, in printCtfData
    defocusratio = ctfvalue['defocus2']/ctfvalue['defocus1']
ZeroDivisionError: float division

Thanks for getting to this, I'll look around to see if I can help diagnose the problem.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Michael Cianfrocco (@mcianfrocco) Original Date: 2014-05-08T16:31:01Z

In looking into the code, the bug from running ctffind comes from getCtfExtrema within ctftools.py. I'm going to walk through my own thinking behind this part of the code, but it would be worth getting it double checked:

This definition is trying to solve a quadratic equation for (frequency)^2 in the CTF equation:

So that when you have the quadratic equation:

Which means that when solving the quadratic equation, you are dividing by 2*a, which is 0.

So, for this special case of Cs = 0, the CTF equation actually changes, allowing us to solve it differently:

Reorganizing and solving for frequency:

Therefore, when Cs = 0, we don't need to solve a quadratic, and can instead solve this simpler case for the radius. So I propose the changes to the code to be: 

def getCtfExtrema(focus=1.0e-6, mfreq=1.498e-04, cs=2e-2,
                volts=120000, ampconst=0.000, numzeros=3, zerotype="peaks"):
        """
        mfreq - frequency in inverse meters = 1.0/(mpix * numcols)
        """
        if debug is True:
                print "defocus %.2f microns (underfocus is positive)"%(focus*1e6)
                print "Freq %.1e 1/m"%(mfreq)
                print "C_s %.1f mm"%(cs*1e3)
                print "High tension %.1f kV"%(volts*1e-3)
        if focus*1e6 > 15.0 or focus*1e6 < 0.1:
                apDisplay.printWarning("atypical defocus value %.1f microns (underfocus is positive)"
                        %(focus*1e6))
        if cs*1e3 > 7.0 or cs*1e3 < 0.0:
                apDisplay.printWarning("atypical C_s value %.1f mm"%(cs*1e3))
        if mfreq > 1e7 or mfreq < 1e5:
                apDisplay.printWarning("atypical mfreq value %.2e 1/meters"%(mfreq))
        if volts*1e-3 > 400.0 or volts*1e-3 < 60:
                apDisplay.printWarning("atypical high tension value %.1f kiloVolts"%(volts*1e-3))

        wavelength = getTEMLambda(volts)

    if cs > 0:
            a = 0.5*cs*math.pi*wavelength**3
            b = -focus*math.pi*wavelength
            c = -math.asin(ampconst)
            if debug is True:
                    print "quadradtic parameters %.3e, %.3e, %.3e"%(a,b,c)
            #eq: sin^2 (a r^4 + b r^2 + c) = 0
            #==> a r^4 + b r^2 + c - n*pi/2 = 0
            #quadradtic: r^2 = [ -b +/- sqrt( b^2 - 4*a*(c + n*pi/2)) ] / 2*a
            # note: "-b + sqrt(..)" is always the positive (non-imaginary) root 

            ## after a certain point the peaks switch direction
            #peakswitch = (2.0*math.sqrt(focus/(cs*wavelength**2)))/math.pi + 0.9
            #if debug is True:
            #       print "Peak switch", peakswitch

            distances = []
            for i in range(numzeros):
                    if zerotype.startswith("valley"):
                            innerroot = b**2. - 4. * a * (c + (i+1)*math.pi)        ## just valleys/minima
                    elif zerotype.startswith("peak"):
                            innerroot = b**2. - 4. * a * (c + (i+0.5)*math.pi)      ## just peaks/maxima
                    else:
                            innerroot = b**2. - 4. * a * (c + (i/2.0+0.5)*math.pi)  ## all extrema
                    if innerroot < 0:
                            continue
                    root = math.sqrt(innerroot)
                    radsq1 = (-b + root)/(2*a)
                    radsq2 = (-b - root)/(2*a)
                    if radsq1 > 0 and radsq1 < radsq2:
                            rad1 = math.sqrt(radsq1)
                            pixeldist = rad1/mfreq
                    elif radsq2 > 0 and radsq2 < radsq1:
                            rad2 = math.sqrt(radsq2)
                            pixeldist = rad2/mfreq
                    else:
                            print "ERROR"
                            continue
                    distances.append(pixeldist)
                    if debug is True:
                            print "radius of zero number %d is %d pixels"%(i+1, pixeldist)

    if cs == 0:

                b = focus*math.pi*wavelength
                c = math.asin(ampconst)

        #I don't understand what the zerotype.startswith('vally' or 'peak') is doing, so I didn't
        #include this.  

        rad = math.sqrt(c/b)        
        pixeldist = rad/mfreq

        distances.append(pixeldist)
                        if debug is True:
                                print "radius of zero number %d is %d pixels"%(i+1, pixeldist)

    return numpy.array(distances)

Where I just added a conditional statement for Cs > 0 or Cs = 0, and provide an equation for the radius.

I don't know if this will be helpful, but I thought I'd start the discussion.

thanks, mike

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Anchi Cheng (@anchi2c) Original Date: 2014-05-08T19:37:42Z

I add Neil to the watcher.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Amber Herold (Amber Herold) Original Date: 2014-05-12T14:58:57Z

I've added Bridget, Clint and David to the discussion. I brought this up in our group meeting this morning, and there were some concerns about using 0 cs. Bridget will take a closer look at this.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Bridget Carragher (@bcarr-czi) Original Date: 2014-05-14T02:08:56Z

Hi all, Is this a legitimate thing to do? i.e. using the standard CTF code with Cs=0? What do others who use Cs correctors do? Has Holger published on this? Seems like it needs to be worked out offline and we need to make sure it all is legit before we try and hack it into Appion? But perhaps I am late tho this party and this has already been discussed and worked out? Bridget

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Neil Voss (@vosslab) Original Date: 2014-05-14T09:59:02Z

I do not think there would be a problem with the code when using Cs = 0. It may be slower than necessary (no Cs gives a simple CTF equation that is more manageable), but I would hope Appion would just work as is.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Michael Cianfrocco (@mcianfrocco) Original Date: 2014-05-14T12:12:12Z

As Neil said, I think that there shouldn't be any theoretical issues with estimating the defocus of images taken at Cs = 0. It is worth pointing out that we are using CTFFIND to estimate the defocus of the images using Cs=0 and it outputs reasonable defocal values without crashing or running into any problems.

So, I think the problem lies within the generation of the CTF output image display, not within the estimation of the defocus. This is because the code was written without taking the possibility of Cs=0 into consideration.

Even though the cryo-EM field hasn't ruled yet on whether spherical aberration corrected microscopes are 'worth it', I think there will be a number of other Appion users who will run into similar problems. For instance, the shared Titan Krios at Janelia Farm is Cs corrected, so I'm assuming this will come up again and again.

For our lab here, I've implemented a Cs override option within ctfestimate.py, allowing the user to specify a Cs value at the command line that is used during CTF estimation but not when it is generating the images for display (which is where the program crashes). While not ideal, it will allow our users to use Appion to analyze their huge datasets from Janelia Farm while this issue is worked out.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Bridget Carragher (@bcarr-czi) Original Date: 2014-05-14T12:13:59Z

Oh OK, got it. Neil can you take a look at why the display code is ding a divide by zero?

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Neil Voss (@vosslab) Original Date: 2014-05-14T18:59:15Z

Oh sorry, I came in late and I did not read the whole thread.

The code above does need to address the @zerotype.startswith('valley' or 'peak')@, because sometimes another program will only request the peaks of the CTF equation, sometimes the valley, and sometimes both, so that has to be addressed or it will crash, as it does. I can look into this tomorrow.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Neil Voss (@vosslab) Original Date: 2014-05-14T19:00:46Z

I will probably need a Cs corrected image to test this on. Michael could you upload a test image to this bug report by chance.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Michael Cianfrocco (@mcianfrocco) Original Date: 2014-05-14T20:06:37Z

Here is a micrograph that is Cs corrected. These are specs:

Nominal Mag: 59,000 Pixel size: 1.05 A/pix

Let me know if you need anything else.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Neil Voss (@vosslab) Original Date: 2014-05-14T22:13:15Z

Hi Michael, what is the high tension voltage. Thanks.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Michael Cianfrocco (@mcianfrocco) Original Date: 2014-05-14T22:59:08Z

Oh, right! It is 300 kV.

leginonbot commented 8 months ago

Original Redmine Comment Author Name: Neil Voss (@vosslab) Original Date: 2014-05-16T15:58:19Z

I just committed r18281 and it works now. You were pretty close Michael, but you needed to return all zeros of the CTF not just one.

I was able to run CTFFIND, PhasorCTF, and RefineCTF, but ACE2 crashed. ACE2 did not like the Cs = 0, probably because it uses the same quadratic equation.